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Abstract. This paper presents results from axisymmetric simulations of magneto-rotational stellar core collapse to 
neutron stars in general relativity using the passive field approximation for the magnetic field. These simulations 
are performed using a new general relativistic numerical code specifically designed to study this astrophysical 
scenario. The code is an extension of an existing (and thoroughly tested) hydrodynamics code, which has been 
applied in the recent past to study relativistic rotational core collapse. It is based on the conformally-flat ap- 
proximation of Einstein's field equations and conservative formulations of the magneto-hydrodynamics equations. 
The code has been recently upgraded to incorporate a tabulated, microphysical equation of state and an approxi- 
mate deleptonization scheme. This allows us to perform the most realistic simulations of magneto-rotational core 
collapse to date, which are compared with simulations employing a simplified (hybrid) equation of state, widely 
used in the relativistic core collapse community. Furthermore, state-of-the-art (unmagnetized) initial models from 
stellar evolution are used. In general, stellar evolution models predict weak magnetic fields in the progenitors, 
which justifies our simplification of performing the computations under the approach that we call the passive 
field approximation for the magnetic field. Our results show that for the core collapse models with microphysics 
the saturation of the magnetic field cannot be reached within dynamical time scales by winding up the poloidal 
magnetic field into a toroidal one. We estimate the effect of other amplification mechanisms including the magneto- 
rotational instability (MRI) and several types of dynamos. We conclude that for progenitors with astrophysically 
expected (i.e. weak) magnetic fields, the MRI is the only mechanism that could amplify the magnetic field on 
dynamical time scales. The uncertainties about the strength of the magnetic field at which the MRI saturates 
are discussed. All our microphysical models exhibit post-bounce convective overturn in regions surrounding the 
inner part of the proto-neutron star. Since this has a potential impact on enhancing the MRI, it deserves further 
investigation with more accurate neutrino treatment or alternative microphysical equations of state. 
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1. Introduction 

Understanding the dynamics of the gravitational collapse 
of the core of massive stars leading to supernova explo- 
sions still remains one of the primary problems in general 
relativistic astrophysics, despite the continuous theoretical 
efforts during the last few decades. This problem stands 
as a distinctive example of a research field where essential 
progress has been accomplished through numerical mod- 
elling with increasing levels of complexity in the input 
physics: hydrodynamics, gravity, magnetic fields, nuclear 
physics, equation of state (EOS), neutrino transport, etc. 
While studies based upon Newtonian physics are highly 
developed nowadays, state-of-the-art simulations still fail, 



broadly speaking, to generate successful supernova explo- 
sions under generic conditions (see e.g. Buras et al. 2003; 
Kifonidis et al. 2006 for details on the degree of sophisti- 
cation achieved in present-day supernova modelling, and 
Woosley & Janka 2005 and references therein for a review 
on the mechanism of core collapse supernovae). The rea- 
sons behind those apparent failures are diverse, all having 
to do with the limited knowledge of some of the underlying 
key issues such as realistic precollapse stellar models (in- 
cluding rotation, or the strength and distribution of mag- 
netic fields), the appropriate EOS, as well as numerical 
limitations due to the need for Boltzmann neutrino trans- 
port, multi-dimensional hydrodynamics, and relativistic 
gravity. 
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Aside from their assistance to understand the super- 
nova mechanism, numerical simulations of stellar core col- 
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lapse are highly motivated nowadays by the prospects of 
a direct detection of the gravitational waves emitted in 
this scenario. In core collapse events where rotation plays 
a role, one of the emission mechanism for gravitational 
waves is the hydrodynamic core bounce, which generates 
a burst signal. The post-bounce wave signal also exhibits 
large amplitude oscillations associated with pulsations in 
the collapsed core (Zwerger & Miiller 1997; Rampp et al. 
1998), neutrino-driven convection behind the supernova 
shock (Miiller et al. 2004) and (possibly) rotational dy- 
namical instabilities (Ott et al. 2005; Ott et al. 2007a,b). 
However, a successful future detection of gravitational ra- 
diation from stellar core collapse faces the combined prob- 
lem of the smallness of the signal strength and of the 
complexity of the burst signal from bounce. On the other 
hand, the energy released in gravitational waves is so small 
that its backreaction to the collapse dynamics is negligible, 
which can significantly simplify the numerical simulation 
of this scenario. To pave the road for a successful detec- 
tion through waveform templates for data analysis, such 
simulations are essential. 

At birth neutron stars have intense magnetic fields 
(<~ 10 12 -10 13 G) or in extreme cases even larger ones 
(~ 10 14 -10 15 G), as inferred from studies of anomalous 
X-ray pulsars and soft gamma-ray repeaters (Kouveliotou 
et al. 1998). For magnetars, the magnetic field can be so 
strong as to alter the internal structure of the neutron 
star (Bocquet et al. 1995). The emergence of such strong 
magnetic fields in neutron stars from the initial field con- 
figuration in the pre-collapse stellar cores is an active and 
important field of research. Similarly, the rotation state 
of the nascent proto-neutron star (PNS) is determined 
by the amount and distribution of angular momentum 
in the core of the progenitor, which is still rather un- 
constrained, being only currently incorporated into stel- 
lar evolution codes (Heger et al. 2005). Observations of 
surface velocities imply that a large fraction of progenitor 
cores is rapidly rotating. The presence of intense mag- 
netic fields, on the other hand, may also affect rotation in 
the core, as it can be spun down in the red giant phase 
by magnetic torques via dynamo action which couples to 
the outer layers of the star (Meier et al. 1976; Spruit & 
Phinney 1998; Spruit 2002; Heger et al. 2005). The lat- 
est numerical calculations of stellar evolution thus predict 
low pre-collapse core rotation rates, which are in agree- 
ment with observed periods of young neutron stars in the 
range of ~ 10-15 ms. Nevertheless, a recent estimate by 
Woosley & Heger (2006) indicates that ~ 1% of all stars 
with M > 10M Q could still have rapidly rotating cores, 
which could also be relevant for the collapsar-type gamma- 
ray burst scenario. 

The presence of intense magnetic fields in a PNS makes 
magneto-rotational core collapse simulations mandatory. 
The weakest point of all existing simulations to date is 
the fact that both the strength and distribution of the 
initial magnetic field in the core are basically unknown. 
If the magnetic field is initially weak, there exist several 
mechanisms which may amplify it to values which can 



have an impact on the dynamics, among them differential 
rotation (fi-dynamo 1 , the magneto-rotational instability 
(MRI hereafter), as well as dynamo mechanisms related 
to convection or turbulence. The first of these mechanisms 
transforms rotational energy into magnetic energy, wind- 
ing up any seed poloidal field into a toroidal field. The 
MRI leads to a sustained turbulent dynamo which is able 
to transport angular momentum outwards, although it re- 
mains unclear how large the actual amplification by this 
process can be (see below). The latter mechanisms, which 
are generically called a-fi-dynamo and will be discussed 
below, include a number of processes which can also pro- 
duce an amplification of the magnetic field. 

Magneto-rotational core collapse simulations were first 
performed as early as in the 1970s (LeBlanc & Wilson 
1970; Bisnovatyi-Kogan et al. 1976; Meier et al. 1976; 
Miiller & Hillebrandt 1979; Ohnishi 1983; Symbalisty 
1984), in which magneto-rotational core collapse was al- 
ready proposed as a plausible supernova explosion mech- 
anism. In recent years, an increasing number of authors 
have performed axisymmetric magneto-hydrodynamic 
(MHD) simulations of stellar core collapse (within the 
so-called ideal MHD limit) employing a Newtonian treat- 
ment of MHD and gravity, and either a simplified equa- 
tion of state (Yamada & Sawai 2004; Ardeljan et al. 2005; 
Sawai et al. 2005) or a microphysical description of matter 
(Kotake et al. 2004a, b, 2005). The main implications of the 
presence of strong magnetic fields in the collapse are the 
redistribution of the angular momentum and the appear- 
ance of jet-like explosions. Specific magneto-rotational ef- 
fects on the gravitational wave signature were first studied 
in detail by Kotake et al. (2004a) and Yamada & Sawai 
(2004), who found differences with purely hydrodynamic 
models only for very strong initial fields (> 10 12 G). The 
most exhaustive parameter study of magneto-rotational 
core collapse to date has been carried out very recently by 
Obergaulingcr et al. (2006a, b). Their axisymmetric sim- 
ulations employed rotating polytropes, Newtonian hydro- 
dynamics and gravity (approximating general relativistic 
effects via an effective relativistic gravitational potential 
in their latter work), and ad-hoc initial poloidal magnetic 
field distributions. These authors found that for weak ini- 
tial fields (< 10 11 G, which is the astrophysically most mo- 
tivated case) there are no differences compared to purely 
hydrodynamic simulations, neither in the collapse dynam- 
ics nor in the resulting gravitational wave signal. However, 
strong initial fields (> 10 12 G) manage to slow down the 
core efficiently (leading even to retrograde rotation in the 
PNS) which causes qualitatively different dynamics and 
gravitational wave signals. For the most strongly magne- 
tized models Obergaulinger et al. (2006b) found highly 
bipolar, jet-like outflows. 

Nowadays, there exists sophisticated numerical tech- 
nology to perform general relativistic hydrodynamics sim- 

1 Note that the "fi-dynamo" is also referred to in the litera- 
ture as "wind-up" or "field- wrapping" . We follow in this paper 
the notation used by Obergaulinger et al. (2006a,b) 
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ulations (see e.g. Font 2003). In recent years this tech- 
nology has been extended to general relativistic magneto- 
hydrodynamics (GRMHD) (Koide et al. 1999; De Villiers 
& Hawley 2003; Del Zanna et al. 2003; Gammie et al. 
2003; Duez et al. 2005; Anton et al. 2006). General rela- 
tivistic simulations involve the challenging computational 
task of solving Einstein's field equations coupled to the 
fluid (and magneto-fluid) evolution. The first general rela- 
tivistic axisymmetric simulations of rotational stellar core 
collapse to neutron stars were performed by Dimmclmeier 
et al. (2001, 2002a, b). These simulations employed simpli- 
fied models to describe the thermodynamics of the pro- 
cess, in the form of a polytropic EOS modified such that 
it accounts both for the stiffening of the matter above 
nuclear density as well as thermal heating by the pass- 
ing shock front (the so-called hybrid EOS; see Janka 
et al. 1993). The inclusion of relativistic effects within 
the so-called CFC approximation results primarily in a 
stronger gravitational pull during the contraction of the 
core. Thus, higher densities than in Newtonian mod- 
els are reached during bounce, and the nascent PNS is 
more compact. Relativistic simulations with improved dy- 
namics and gravitational waveforms were reported by 
Cerda-Duran et al. (2005), who used the CFC+ frame- 
work, which includes second post-Newtonian corrections 
to CFC. Comparisons of the CFC approach with fully 
general relativistic simulations (employing also stable re- 
formulations of the Einstein equations in 3 + 1 form) have 
been reported by Shibata & Sekiguchi (2004), Ott et al. 
(2007a), and Ott et al. (2007b) in the context of axisym- 
metric core collapse simulations. As in the case of CFC+, 
the differences found in both the collapse dynamics and 
the gravitational waveforms arc minute, which highlights 
the suitability of CFC for performing accurate simula- 
tions of such scenarios without the need for solving the 
full system of Einstein's equations. Owing to the excellent 
approximation of full general relativity offered by CFC 
in the context of stellar core collapse, extensions to im- 
prove the microphysics through the incorporation of a tab- 
ulated non-zero temperature EOS and a simplified neu- 
trino treatment have been recently reported by Ott et al. 
(2007a) and Dimmclmeier et al. (2007a). These simula- 
tions allow a direct comparison to the models presented 
in Dimmclmeier et al. (2002b), Cerda-Duran et al. (2005), 
and Shibata & Sekiguchi (2004), which use the same pa- 
rameterization of rotation but a simple hybrid EOS. This 
comparison shows that with a microphysical treatment the 
influence of rotation on the collapse dynamics and wave- 
forms is significantly reduced. In particular, the most im- 
portant result of Dimmclmeier et al. (2007a) is the sup- 
pression of core collapse with multiple centrifugal bounces 
and its associated Type II gravitational waveforms (see 
Dimmclmeier et al. 2002b). 

On the other hand, to further improve the realism of 
core collapse simulations in general relativity, the incor- 
poration of magnetic fields in numerical codes via solving 
the MHD equations is also currently being undertaken 
(Shibata et al. 2006; Cerda-Duran & Font 2007). The 



work of Shibata et al. (2006) is focused on the collapse 
of initially strongly magnetized cores (~ 10 12 -10 13 G). 
Although these values are probably astrophysically not 
relevant (as the stellar evolution models of Heger et al. 
2005, predict a poloidal field strength of <~ 10 6 G in the 
progenitor), they enable them to resolve the scales needed 
for the MRI to develop, since the MRI length scale grows 
with the magnetic field. The results of Shibata et al. (2006) 
show that the magnetic field is mainly amplified by the 
wind-up of the magnetic field lines by differential rotation. 
Consequently, the magnetic field is accumulated around 
the inner region of the PNS, and a MHD outflow forms 
along the rotation axis removing angular momentum from 
the PNS. A different approach is followed by Cerda-Duran 
& Font (2007) . Their progenitors are chosen to be weakly 
magnetized (< 10 G) which is in much better agreement 
with predictions from stellar evolution. Under these condi- 
tions the so-called "passive" magnetic field approximation 
(see Sect. 2.3 below) is appropriate. In addition, the nu- 
merical code used in that work, which utilizes spherical 
coordinates, is more suitable for core collapse simulations 
than codes based on Cartesian/cylindrical coordinates, as 
used e.g. by Shibata et al. (2006). 

In this paper we continue the program initiated in 
Cerda-Duran & Font (2007) to build a numerical code 
which includes all relevant ingredients to study relativis- 
tic magneto-rotational stellar core collapse. To this aim we 
present here the first relativistic simulations of magneto- 
rotational core collapse which take into account the effects 
of a microphysical EOS and a simplified neutrino treat- 
ment. Those effects have been incorporated in the code 
following the approach recently presented by Ott et al. 
(2007a) and Dimmclmeier et al. (2007a). As in Cerda- 
Duran & Font (2007) we employ the passive magnetic field 
approximation in the treatment of the magnetic field. 

The paper is organized as follows: Sect. 2 presents a 
brief overview of the theoretical framework we use to per- 
form relativistic simulations of core collapse. Sect. 3 de- 
scribes how the magnetized initial models for core collapse 
are built and presents our sample of models. In Sect. 4 
we discuss aspects related to incorporating microphysics 
in the core collapse models and their implementation in 
the numerical code. A brief outline of our numerical ap- 
proach is discussed in Sect. 5. The evolution of the core 
collapse initial models is discussed in Sect. 6. The main 
paper closes with a summary in Sect. 7. Relevant tests of 
the code are analyzed in Appendix A, while Appendix B 
provides an estimate for the growth rate of the f2-dynamo. 

Throughout the paper we use a spacelike signature 
(—,+,+,+) and units in which c = G = 1. Greek in- 
dices run from to 3, Latin indices from 1 to 3, and we 
adopt the standard Einstein summation convention. 

2. Theoretical framework 

We adopt the 3 + 1 formalism of general relativity 
(Lichnerowicz 1944) to foliate the spacetime into space- 



4 



Cerda-Duran et al.: Relativistic passive-magneto-rotational core collapse 



like hypcrsurfaces. In this approach the line element reads 

d.s 2 = -a 2 dt 2 + (dx l + L At) (dx j + (3 3 dt), (1) 

where a is the lapse function, fi l is the shift vector, and 7^ 
is the spatial three-metric induced in each hypersurface. 
Using the projection operator _L£J and the unit four-vector 
ra M normal to each hypersurface, it is possible to build the 
quantities 



2rp00 



E = n^n v T^ = a z T 

C I M I V rp _ rp 



(2) 
(3) 
(4) 



which represent the total energy, the momenta, and the 
spatial components of the stress-energy tensor, respec- 
tively. 

To solve the gravitational field equations we choose the 
ADM gauge in which the three-metric can be decomposed 
as = </> 4 7ij + hJ^ T , where <f> is the conformal factor, 
7ij is the fiat three-metric, and hJ 3 T is the transverse and 
traceless part of the three-metric. Note that this gauge 
choice implies the maximal slicing condition in which the 
trace K of the extrinsic curvature tensor Ka vanishes. 



2.1. The CFC approximation 

In our work Einstein's field equations are formulated and 
solved using the conformally flat condition (CFC here- 
after), introduced by Isenberg (1978) and first used in a 
dynamical context by Wilson et al. (1996). In this approx- 
imation, the three-metric in the ADM gauge is assumed 
to be conformally flat, 7^ = (f) 4 ^. Note that the same 
aproximation can be realized for other gauge choices such 
as the quasi-isotropic gauge or the Dirac gauge, both sup- 
plemented by the maximal slicing condition. Under the 
CFC assumption the gravitational field equations can be 
written as a system of five nonlinear elliptic equations, 



1D7T 



7K-K l i 

AM) = 27ra0 5 E + 2S + 10 

16-7T 



(5) 
(6) 



A/3* = lQirac^S 1 



^K^j -lyv k {3\ (7) 



where A and V are the Laplace and nabla operators as- 
sociated with the flat three- metric, and S — 7^ Sij . 



the magnetic field b^ measured by a comoving observer. 
In this case the total energy-momentum tensor is given by 



= (ph + b 2 ) U»U» + (P+—)g^- b^b", 



(8) 



where p is the rest- mass density, h = l + e + P/p the rela- 
tivistic enthalpy, e the specific internal energy, P the pres- 
sure, and the four-velocity of the fluid, while b 2 = b^b^. 
We define the magnetic pressure P mag = b 2 /2 and the spe- 
cific magnetic energy e mag = b 2 /(2p), whose effect on the 
dynamics is similar to that played by the pressure and spe- 
cific internal energy of the fluid, respectively In the ideal 
MHD limit, the electric field measured by a comoving ob- 
server vanishes, and Maxwell's equations simplify. Under 
this assumption the electric field four-vector can be 
expressed in terms of the magnetic field four- vector _B M , 
and only equations for B % are needed. For an Eulerian ob- 



server, 



the temporal component of the electric 



field vanishes, E 1 * = (0, —EijkV 3 B k ). In this case Maxwell's 
equations reduce to the divergence-free condition and the 
induction equation for the magnetic field, 



ViB* 1 = 0, 



= W j (v* i B* j - v* ] B*' 



(9) 



where B* % — \/jB l and v*' 1 — av l — f3 l , with v l being 
the fluid's three-velocity as measured by the Eulerian ob- 
server. The ratio of the determinants of the three-metric 
and the flat three-metric is given by 7 = 7/7. In the 
Newtonian limit v* 1 — ► v 1 and B* 1 — > B\ and the 
Newtonian induction equation and divergence constraint 
are recovered. 

The evolution of a magnetized fluid is determined by 
the conservation law of the energy-momentum, V^T 7 ^ = 
0, and by the continuity equation, V^J^ 1 = 0, for the 
rest-mass current J M = pu^ . Following the procedure laid 
out in Anton et al. (2006), in order to cast the GRMHD 
equations as a hyperbolic system of conservation laws well 
adapted to numerical work, the conserved quantities are 
chosen in a way similar to the purely hydrodynamic case 
presented by Banyuls et al. (1997): 



D = pW, 

5, = (ph + b 2 )W 2 v l -ahb , 



T = (ph + b 2 )W 2 



0\2 



a z (b K> ) 



D, 



(10) 
(11) 
(12) 



2.2. General relativistic magnetohydrodynamics 

The energy-momentum tensor of a magnetized perfect 
fluid can be written as the sum of the fluid part and 
the electromagnetic field part. In the so-called ideal MHD 
limit (where the fluid is a perfect conductor of infinite con- 
ductivity), the latter can be expressed solely in terms of 



where W = au° is the Lorentz factor. With this choice 
the system of conservation equations for the fluid and the 
induction equation for the magnetic field can be cast as a 
first-order, flux-conservative, hyperbolic system, 



dt 



dx % 



(13) 
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with the state vector, flux vector, and source vector given 

by 

U={D,Sj,r,B k ], 



W 



Dv\S j V + 5)[P+- 
tv 1 + (p + ^ v i - a^-,v t B k - v k B l 



0, -T M ' 

' 2 dxi 



flV J 1 



(14) 

(15) 
(16) 



where <5* is the Kronecker delta and are the Christoffcl 
symbols associated with the four-metric. We note that the 
above definitions contain components of the magnetic field 
measured by both a comoving observer and an Eulerian 
observer. The two are related by 



WB l Vi 



a 



b l = 



B l + ab°i 
W 



The hyperbolic structure of Eq. (13) and the associated 
spectral decomposition (into eigenvalues and eigenvectors) 
of the flux- vector Jacobians is given in Anton et al. (2006). 
This information is needed for numerically solving the sys- 
tem of equations using the class of high-resolution shock- 
capturing schemes that we have implemented in our code. 

2.3. The passive field approximation 

In the collapse of weakly magnetized stellar cores, it is a 
fair approximation to assume that the magnetic field en- 
tering in the energy-momentum tensor of Eq. (8) is negli- 
gible when compared with the fluid part, i.e. P mag <C P, 
Cmag "C e, and that the components of the anisotropic term 
of T» v satisfy WW < phuV + Pg^ v . With these simpli- 
fications the evolution of the magnetic field, governed by 
the induction equation, does not affect the dynamics of the 
fluid, which is governed solely by the hydrodynamics equa- 
tions. However, the magnetic field evolution does depend 
on the fluid evolution, due to the presence of the velocity 
components in the induction equation. This "test mag- 
netic field" (or passive field) approximation is employed 
in the core collapse simulations reported in this work. 

Within this approach the seven eigenvalues of the 
GRMHD Riemann problem (entropy, Alfven, and fast and 
slow magnetosonic waves) reduce to three (Cerda-Duran 
& Font 2007), 



A. 



hydro 



A! = A 



A± 



A, 



A ± hydro ~ A f ± i 



(18) 
(19) 



where AQ hydro and A! thydro are the eigenvalues of the 
Jacobian matrices of the purely hydrodynamics equations, 
as reported by Banyuls et al. (1997). 

This approximation has several interesting properties. 
First, if we perform a simulation for a given initial mag- 
netic field, we can compute the result for a simulation 
with the same initial magnetic field scaled by some fac- 
tor. To do this it is sufficient to increase or reduce the 



strength of the magnetic field at any given time during 
the simulation by the same factor. The second property 
is what we call the "composition rule" . If we perform two 
simulations with the same hydrodynamics but different 
initial magnetic fields, B\ Q and B^ , any linear combina- 
tion B*(t) = a B\(t) + b B^if) of the magnetic field at any 
time, with a and b being constants, will be the solution for 
the evolution of a model whose initial magnetic field is the 
same linear combination, i.e. B* a = aB^+bB^ . Hence, 
we can make use of these properties to cover a wide range 
of magnetic field strengths and structures by performing 
just a few simulations, and then constructing additional 
ones by means of the "composition rule" . Needless to say, 
these two properties are valid only if the magnetic field re- 
sulting from the scaling or composition satisfies itself the 
passive field approximation at all times. 



(17) 2.4. Gravitational waves 



The Newtonian standard quadrupole formula has been ex- 
tensively used in numerical simulations of astrophysical 
systems to compute the gravitational radiation and wave- 
forms without having to consider the full evolution of the 
spacetime and solving Einstein's equations. This formula 
computes the radiative part of the spatial metric as 



quad 



P 



TTfci 



R ~* ' 



(20) 



where P^j™ is the transverse traceless projector operator 
(Thorne 1980), R is the distance to the source, Qij is 
the mass quadrupole moment, and a dot denotes a time 
derivative. In spite of its simplicity, the particular form in 
which Eq. (20) is expressed leads to numerical difficulties 
due to the presence of second time derivatives. A way to 
circumvent this problem is to eliminate all time derivatives 
using the equations of motion. Following Finn (1989) and 
Blanchet et al. (1990) one can arrive to an expression for 
Qki with no explicit appearance of time derivatives. This 
is the so-called stress formula, 

Qij « STF 1 2 J d 3 x ^D* v k v l + x k j kl V, C/) J , 

(21) 

where STF means the symmetric and traceless part, and 
D* = V7-D. This formula has proved to be numerically 
much more accurate than the original formula and we use 
it in this paper to extract gravitational waveforms. 

In the case of a magnetized fluid in the ideal MHD case, 
the gravitational radiation is also affected by the energetic 
content of the magnetic field. Kotake et al. (2004b) have 
derived an extension of the quadrupole formula for such 
a case. In a similar way, it is possible to calculate the 
corresponding stress formula (Obergaulinger et al. 2006b), 
which reads 



Qij « STF<^ 2 / d 6 x V7 



D* (liklji v k v l + x k % t VjU} 



-likljl b k b l 



(22) 
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Note that in the limit of weak magnetic fields the orig- 
inal stress formula is recovered. We use this formula in 
the magnetized core collapse simulations to calculate the 
contribution of the magnetic field to the waveforms in the 
passive field approximation. 

3. Initial data 

The structure and strength of the magnetic field in the 
stellar core collapse progenitors, needed as initial condi- 
tions of our numerical simulations, is still an open ques- 
tion in astrophysics. State-of-the-art models from stellar 
evolution including a description for the influence of the 
magnetic field (Heger et al. 2005), predict that the distri- 
bution of the magnetic field in the iron core has proba- 
bly a dominant toroidal component, with a strength of 
about 10 9 -10 10 G, and a poloidal component of only 
about 10 5 -10 6 G. For such weak fields (P mag < P), 
the passive field approximation adopted here is likely to 
be sufficient to describe the initial models considered in 
this work. A second consideration is whether or not the 
initial model should be an equilibrium model. In gen- 
eral, if one tries to construct a stationary model without 
meridional currents and assuming an iscntropic flow, the 
only possible magnetic field configuration is poloidal (see 
Bekenstein & Oron 1979). Stationary models of magne- 
tized stars have been computed under these assumptions 
by Bocquet et al. (1995). In the general (but still isen- 
tropic) case in which meridional circulation is allowed, a 
toroidal component of the magnetic field may exist, but 
the method to calculate stationary models is far more com- 
plicated (Gourgoulhon & Bonazzola 1993; Ioka & Sasaki 
2003, 2004). When one considers ideal MHD, also purely 
toroidal magnetic fields exist which maintain the circu- 
larity condition (Oron 2002), and therefore it is possible 
to generate stationary models without meridional compo- 
nents. Finally, in the case that magnetic pressure does not 
exceed the hydrostatic pressure, Oron (2002) has shown 
that stationary models with mixed toroidal and poloidal 
component approximately accomplish the circularity con- 
dition. 

Therefore, it makes sense to construct initial models 
for magnetized stellar cores by simply adding an ad-hoc 
magnetic field to a purely hydrodynamic equilibrium con- 
figuration. If the condition B* ■ VfT = is satisfied, where 
O* = v* v /(r sin 6) is the angular velocity of the fluid, then 
the initial magnetic field does not evolve in time either. 
Note that in this work we use as initial models both equi- 
librium and non-equilibrium configurations for the mag- 
netic field, as specified in Table 2. 

3.1. Magnetic field configurations 

Since the numerical scheme we use to evolve the MHD 
equations only preserves the value of V • B* but does 
not impose the divergence constraint of the magnetic field 
itself, it is necessary to build initial configurations that 
also satisfy this condition. To do this we calculate the 



initial magnetic field from a vector potential A* , such that 
B* = V x A*, which is discretized as explained in Cerda- 
Duran & Font (2007). 

For our code tests and core collapse simulations we 
use three possible magnetic field configurations as initial 
conditions (or any possible combination of them): 

- the homogeneous "starred" magnetic field, in which B* 
is constant and parallel to the symmetry axis, 

- the poloidal magnetic field generated by a circular cur- 
rent loop of radius r mag (Jackson 1962), that can be calcu- 
lated from the only non- vanishing component of the vector 
potential A* as 

A <P ~ 2 ^ 2 n (n + lV 2n+2^2n+lK COS V)' 

(23) 

where r< = min(r, r mag ), r> = max(r, r mag ), and Bq is 
the magnetic field at the center, and 

- a toroidal magnetic field of the form 

r 2 

B**=B* mag (24) 

' mag M 

whose maximum value is reached at w = r mag , where Bq 
is the initial central magnetic field and w is the distance 
to the axis. 

Note that in all three cases, we employ the "starred" 
magnetic field, since the divergence constraint is valid for 
this quantity when computed with respect to the flat di- 
vergence operator. In this way we can extend any analytic 
prescription for the magnetic field given in flat spacetime 
in an easy way. Also note that in the presence of strong 
gravitational fields the magnetic field B is deformed with 
respect to B* due to the curvature of the spacetime, al- 
though the divergence constraint is automatically fulfilled. 
Fig. 1 shows examples of the magnetic field structure for 
the poloidal configurations of the initial magnetic field. 

The initial magnetic field configuration is denoted in 
the names of the models in our sample by adding a label to 
the purely hydrodynamic model. For the latter we follow 
the notation of Dimmelmeier et al. (2002a) and Ott et al. 
(2007a). These models are listed in Table 1. The label 
for the magnetic field is constructed following the nota- 
tion of Obergaulinger et al. (2006b). We add the suffix 
D3M0 to denote those models with purely poloidal mag- 
netic field generated by a circular loop and r mag = 400 km 
(M0 denotes the passive field approximation) and we use 
the suffix T3M0 for models with purely toroidal magnetic 
field and r mag = 400 km. We have also built the DT3M0 
model, whose magnetic field distribution is a combination 
of D3M0 and T3M0 with equal magnetic field strengths 
at the center. This model allows us to check the validity 
of the "composition rule" (see Sect. 2.3). 
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Fig. 1. Logarithm of the rest mass density (log p, color coded) and magnetic field lines (B l , white curves) in a homogeneous B* 1 
field configuration (left panel) and a magnetic field generated by a circular current loop (right panel) , for a typical rotating iron 
core (model A1B3 in Table 1) used for the collapse simulations. The magnetic field lines have been calculated as isocontours of 
the vector potential. The axis scale is in km and the density in cgs units. 



3.2. Hydrodynamic configurations 

3.2.1. Polytropes in rotational equilibrium 

For the simulation with a simplified description of matter 
using the hybrid EOS (see Sect. 4.1), we construct 7i„i = 
4/3 polytropes in rotational equilibrium which we obtain 
by using the relativistic generalization of Hachisu's self- 
consistent field method by Komatsu et al. (1989) 2 . Their 
rotation law for the specific angular momentum j is given 

by 



j = A 2 (fl c - fi), 



(25) 



where A parametrizes the degree of differential rotation 
(stronger differentiality with decreasing A) and fi c is the 
value of the angular velocity O at the center. In the 
Newtonian limit, this reduces to 



A 2 Q r 



A 2 



(26) 



The parameters of the selected models, which are cho- 
sen to be identical to some of the models considered by 
Dimmelmeier et al. (2002a), are described in Table 1. In 
addition, as we aim at comparing our results with the re- 
cent numerical simulations performed by Obergaulinger 
et al. (2006b) in Newtonian gravity, a subset of our mod- 
els (those with purely poloidal magnetic field) have been 
selected as general relativistic counterparts of their mod- 
els. In Table 1 we also give the values for the gravitational 
mass M g (which is identical to the ADM mass Madm) and 
for the initial rotation rate (3 — E Toi /\Eb\. In the defini- 
tion of (3 we use the following expressions for the rotational 



The adiabatic index should not be confused with the de- 
terminant of the spacetime three-metric, although we use the 
same symbol 7 (following usual practice). 



Table 1. Purely hydrodynamic initial models used in the mag- 
netized core collapse simulations. 



Model 




A 


P 


M g 


7i 




[f0 10 gem" 3 ] 


[10 3 km] 


[%] 


[M ] 




AfB3G3 


1.00 


50.0 


0.90 


1.46 


1.300 


AfB3G5 


1.00 


50.0 


0.90 


1.46 


1.280 


A3B3G5 


1.00 


0.5 


0.90 


1.46 


1.280 


A2B4G1 


1.00 


1.0 


1.80 


1.50 


1.325 


A4B5G5 


1.00 


0.5 


4.00 


1.61 


1.280 


S20A1B1 


0.88 


50.0 


0.25 


1.58 




S20A1B5 


0.88 


50.0 


4.00 


1.58 




s20A2B2 


0.88 


1.0 


0.50 


1.58 




s20A2B4 


0.88 


1.0 


1.80 


1.58 




s20A3B3 


0.88 


0.5 


0.90 


1.58 




E20A 


0.42 




0.37 


2.00 





kinetic energy E Iot , the gravitational binding energy E^, 
and the magnetic energy -E mag : 



M g - A/ p - E Iot - E, 



E ro t — 

E h = M. 



mag ■ 



(27) 
(28) 
(29) 



where M p is the proper mass. 

For these simplified initial models the gravitational col- 
lapse is initiated by slightly decreasing the adiabatic index 
from its initial value to 71 < 7; n i = 4/3, which results in 
a loss of pressure support. If no pressure reduction were 
imposed, the purely hydrodynamic initial models would 
remain stationary. However, even in that case the associ- 
ated initial magnetic field may not remain stationary (see 
also Table 2 below). Only a purely toroidal initial mag- 
netic field would not evolve in time, while any magnetic 
field configuration of initial models labeled Al would still 
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stay approximately stationary, since these models rotate 
almost rigidly, and thus the initial magnetic field cannot 
wind up itself strongly. 

3.2.2. Presupernova models from stellar evolution 
calculations 

As initial models for the simulations where we use a 
microphysical EOS and deleptonization, we employ the 
solar-metallicity 20 M Q (zero-age main sequence) model 
of Wooslcy et al. (2002) (labeled as model s20 in Tabic 1). 
On this spherically symmetric model, which is initially not 
in equilibrium as it has a non-zero radial velocity profile, 
we impose the rotation law (25), using the same rotation 
nomenclature as for the previously described polytropes 
in equilibrium. 

In addition, we perform calculations with the "rotat- 
ing" presupernova model E20A of Heger et al. (2000), 
which we map onto our computational grids under the 
assumption of constant rotation on cylindrical shells of 
constant vj. 

4. Treatment of matter during the evolution 

In this work we improve upon preceding relativistic stel- 
lar core collapse simulations by using an advanced descrip- 
tion of microphysics as presented in Ott et al. (2007a) and 
Dimmclmeicr et al. (2007a). For comparison to previous 
results, we also perform simulations with the simplified 
hybrid EOS (Janka et al. 1993). In the following, we de- 
scribe both approaches for the treatment of matter. 

4.1. Hybrid EOS 

For calculations employing polytropes in rotational equi- 
librium as initial models, we utilize the hybrid EOS. Here 
the pressure consists of a polytropic part, P p = Kp 1 , with 
K = 4.897 x 10 14 (in cgs units), plus a thermal part, 
Pth — /Q£th(7th — 1), where e t h = e — e P and where we set 
7th = 1.5. The thermal contribution is chosen to take into 
account the rise of thermal energy due to shock heating. 
As p reaches nuclear density at p nuc = 2. Ox 10 14 g cm -3 , 7 
is raised to 72 = 2.5 and K adjusted accordingly to main- 
tain monotonicity of P and e. Due to this stiffening of the 
EOS the core undergoes a so-called pressure-supported 
bounce. More details of the hybrid EOS can be found e.g. 
in Dimmelmeier et al. (2002a). 

4.2. Microphysical EOS, deleptonization scheme, and 
neutrino pressure 

In our more realistic calculations, for which the models 
s20 and E20A from stellar evolution are taken as initial 
models, we employ the tabulated non-zero temperature 
nuclear EOS by Shen et al. (1998) in the variant of Marek 
et al. (2005) which includes baryonic, electronic, and pho- 
tonic pressure components. This gives the fluid pressure 



P (and additional thermodynamic quantities) as a func- 
tion of p, the temperature T, and the electron fraction Y e . 
Since the code operates with the specific internal energy 
e instead of the temperature T, we determine the corre- 
sponding value for T iteratively with a Newton-Raphson 
scheme. 

To determine the evolution of Y e , the state vector, 
flux vector, and source vector for the conservation equa- 
tions (13), as given in Eqs. (14-16) have to be augmented 
by the components 

DY e , DY e v\ S Yc , (30) 

respectively. The source term Sy e is a consequence of the 
electron captures during collapse, which reduces Y e . This 
deleptonization also effectively decreases the size of the 
homologously collapsing inner core, and has thus a direct 
influence on the collapse dynamics and the gravitational 
wave signal. Hence, it is essential to include (at least an 
approximate scheme for) deleptonization during collapse. 

Since multi-dimensional radiation hydrodynamics cal- 
culations in general relativity are not yet computationally 
feasible, in the simulations using the microphysical EOS 
we make use of a a recently proposed scheme (Liebendorfer 
2005) where deleptonization is parametrized based on 
data from detailed spherically symmetric calculations with 
Boltzmann neutrino transport. As in Dimmelmeier et al. 
(2007a) we take the latest available electron capture rates 
(Langankc & Martmcz-Pincdo 2000), which result in lower 
values for Y e in the inner core at bounce compared to 
recent results (Ott et al. 2007a) where standard cap- 
ture rates were used (Rampp & Janka 2000). Following 
Liebendorfer (2005), deleptonization is stopped at core 
bounce (i.e. as soon as the specific entropy s per baryon 
exceeds 3fce). After core bounce Y e is only passively 
advected, neglecting any further deleptonization in the 
nascent PNS. 

Neutrino pressure is included only in the regime which 
is optically thick to neutrinos, which we define for p be- 
ing above the trapping density p t = 2 x 10 12 g cm~ 3 . 
Following Liebendorfer (2005), here we approximate the 
contribution of the neutrino pressure P„ as an ideal Fermi 
gas and include the radiation stress via additional source 
terms in the momentum and energy equations for the fluid. 

5. Outline of the numerical approach 

The GRMHD numerical code we use in our simulations 
is based on the purely hydrodynamic code described in 
Dimmelmeier et al. (2002a, b), and on its extension dis- 
cussed in Cerda-Duran et al. (2005). It has been described 
in detail in a previous paper (Cerda-Duran & Font 2007), 
which allows us to provide here only succint information. 
The code performs the coupled time evolution of the equa- 
tions governing the dynamics of the spacetime, the fluid, 
and the magnetic field in general relativity. The equations 
are implemented in the code using spherical polar coor- 
dinates {t,r,9,(p}, assuming axisymmetry with respect 
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to the rotation axis and equatorial plane symmetry at 
9 = tt/2. 



5.1. The hydrodynamics solver 



For the evolution of the matter fields we utilize a high- 
resolution shock-capturing (HRSC) scheme, which numer- 
ically integrates the subset of equations in system (13) 
that corresponds to the purely hydrodynamic variables 
(D, Si, t). HRSC methods ensure the numerical conser- 
vation of physically conserved quantities and a correct 
treatment of discontinuities such as shocks (see e.g. Font 
2003, for a review and references therein). We have im- 
plemented in the code various cell-reconstruction proce- 
dures, either second-order or third-order accurate in space, 
namely minmod, MC, and PHM (see Toro 1999, for def- 
initions). The time update of the state vector U is done 
using the method of lines in combination with a second- 
order accurate Runge-Kutta scheme. The numerical fluxes 
at the cell interfaces are obtained using either the HLL 
single-state solver of Harten et al. (1983) or the symmet- 
ric scheme of Kurganov & Tadmor (2000) (KT hereafter). 
Both solvers yield results with an accuracy comparable 
to complete Riemann solvers (with the full characteristic 
information), as shown in simulations involving purely hy- 
drodynamic special relativistic flows (Lucas-Serrano et al. 
2004) and general relativistic flows in dynamical space- 
times (Shibata & Font 2005). Tests of both solvers in 
GRMHD have been reported recently by Anton et al. 
(2006). 



5.2. Evolution of the magnetic field 



The evolution of the magnetic field needs to be performed 
in a way that is different from the rest of the conservation 
equations, since the physical meaning of the correspond- 
ing conservation equation is different. Although the induc- 
tion equation can be written in a flux conservative way, 
a supplementary condition for the magnetic field has to 
be given (the divergence constraint), which has to be ful- 
filled at each time iteration. The physical meaning of these 
two equations is the conservation of the magnetic flux in a 
close volume, in our case each numerical cell. Therefore, an 
appropriate numerical scheme has to be used which takes 
full profit of such a conservation law. Among the numeri- 
cal schemes that satisfy this property (see Toth 2000, for 
a review), the constrained transport (CT) scheme (Evans 
& Hawley 1988) has proved to be adequate to perform 
accurate simulations of magnetized flows. Our particular 
implementation of this scheme (see Cerda-Duran & Font 
2007, for details) has been adapted to the spherical polar 
coordinates used in the code. The discretized evolution 
equations for the poloidal components of the magnetic 



field read 



d t B*Zi ■ 



fsin 6 E, 



- [sm6E*] i+hj _ h 



dtB*Li = 2 



r i+ i j A(cose) j 



Ar 2 



(31) 
(32) 



where (in vectorial form) E* = v* x B* , and where cell 
centers are located at (ij), radial interfaces at (i + 
angular interfaces at {ij + |), and cell corners (cell edges 
along the ^-direction) at (i + \j + |). We note that the 
evolution equation for the toroidal magnetic field is analog 
to that used for the hydrodynamics, since in axisymmetry 
this component does not play any role in the CT scheme. 
The previous expressions are used in the numerical code 
to update the magnetic field. The only remaining aspect 
is to give an explicit expression for the value of E* . A 
practical way to calculate E* v from the numerical fluxes in 
the adjacent interfaces (Balsara & Spicer 1999) is 



i+5 j+1 



(33) 



where the fluxes (15) are obtained in the usual way by 
solving the Riemann problem at the interfaces. The com- 
bination of the CT scheme and this way of computing the 
electric field is called the flux-CT scheme. It is used in all 
numerical simulations reported in this paper. Finally, the 
time discretization of Eqs. (31) and (32) is performed in 
the same way as for the fluid evolution equations. 

5.3. The metric solver 

The CFC metric equations are five nonlinear elliptic cou- 
pled Poisson-like equations which can be written in com- 
pact form as Au(x) = f(x;u(x)), where u = u k = 
(4>, a<p, (CP ), and / = f k is the vector of the respective 
sources. These five scalar equations for each component 
of u couple to each other via the source terms that in 
general depend on the various components of u. We use 
a fix-point iteration scheme in combination with a linear 
Poisson solver to solve these equations. Further details on 
this type of metric solver can be found in Cerda-Duran 
et al. (2005) and Dimmclmeier et al. (2002a). 

5.4. Setup of the numerical grid 

We perform all axisymmetric simulations with a resolu- 
tion (n r x no) of 300 x 30 zones, except for models la- 
beled A4B5G5 in which a resolution of 300 x 40 is used 
due to the more complex angular structure. In both cases 
the radial grid is equally spaced for the first 100 points, 
with a grid spacing of 100 m. The remaining radial zones 
are logarithmically distributed to cover the outer parts of 
the star and the exterior artificial low-density atmosphere. 
The angular grid is equally spaced and we assume equato- 
rial symmetry. We have performed resolution tests and we 
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Table 2. Hydrodynamical and magnetic field properties of all models computed in this work. From left to right the columns 
report the name of the model, the stationarity properties of the initial magnetic field, the maximum rest-mass density p max , 
the maximum poloidal field |-B po io|max, the maximum toroidal field l-B^lmax, the rotational energy parameter f3 TOt , the magnetic 
energy parameter /3 ma g, its poloidal contribution /3 po io, and the central angular velocity Q c . The time scale tq of the growth of 
the magnetic field by the f2-dynamo and the saturation time t sat for this process are also shown. 



Model 


stationary 


pmax 


1 R 1 

|-Dpolo |max 


\-t3tp (max 


a 

Prot 


r> 

Pmag 


ft 

Ppolo 


( \ 


tu 


^sat 






r i4 e i 
[ 10 cms] 


[10 10 G] 


[10 10 G] 


[io- 2 ] 


[io- 8 ] 


[io- 8 ] 


[ms- 1 ] 


[s] 


[s] 


A2B4G1-D3M0 


no 


0.47 


400 


1467 


15.6 


8.3 


0.9 


0.36 


85 


10.3 


A1B3G3-D3M0 


approx. 


4.22 


1719 


2522 


2.3 


1.2 


0.8 


3.96 


7 


0.3 


A1B3G3-T3M0 


yes 


4.22 





1714 


2.3 


0.2 


0.0 


3.96 






A1B3G5-D3M0 


approx. 


4.57 


1146 


1275 


0.9 


0.5 


0.4 


3.91 


21 


0.7 


A1B3G5-T3M0 


ves 


4.57 





1542 


0.9 


0.2 


0.0 


3.91 






A1B3G5-DT3M0 


approx. 


4.57 


1146 


1537 


0.9 


0.6 


0.4 


3.91 


21 


0.6 


A3B3G5-D3M0 


no 


3.73 


984 


1672 


2.3 


0.6 


0.4 


3.75 


24 


1.1 


A4B5G5-D3M0 


no 


1.74 


1094 


2716 


8.5 


4.4 


1.2 


1.18 


24 


2.1 


A4B5G5-T3M0 


yes 


1.74 





1626 


8.5 


0.4 


0.0 


1.18 






S20A1B1-D3M0 


approx. 


2.69 


1221 


162 


0.6 


2.9 


2.9 


1.34 


22 


0.5 


s20A2B2-D3M0 


no 


2.75 


1849 


3574 


5.8 


7.6 


3.2 


3.55 


7 


0.5 


s20A2B2-T3M0 


yes 


2.75 





1365 


5.8 


0.9 


0.0 


3.55 






S20A1B5-D3M0 


approx. 


2.69 


1100 


1011 


7.8 


3.2 


2.4 


3.89 


10 


0.8 


S20A1B5-T3M0 


yes 


2.69 





1447 


7.8 


1.2 


0.0 


3.89 






E20A-D3M0 


no 


2.29 


2343 


7503 


7.7 


23.4 


1.8 


4.37 


6 


0.5 


E20A-T3M0 


yes 


2.29 





2739 


7.7 


1.9 


0.0 


4.37 







have found that such a resolution is adequate for our simu- 
lations (see Cerda-Duran 2006; Cerda-Duran & Font 2007, 
for details). As a consequence of our various code tests (see 
Appendix A) all results discussed in Sect. 6 correspond to 
simulations performed using PHM reconstruction and the 
HLL solver for the hydrodynamics. 

6. Results 

We now present the main results from our numerical sim- 
ulations of rotational magnetized core collapse to neutron 
stars. First, we note that a quantitative summary of our 
findings is reported in Table 2, to which we will refer 
repeteadly. The dynamics of the models we have selected is 
identical to the dynamics of the unmagnetized ones, since 
the passive field approximation is used. Therefore, we will 
not describe here all the morphological features of the hy- 
drodynamics of both models with the hybrid EOS (sim- 
plified models hereafter) and models with the microphys- 
ical EOS and the delcptonization scheme (microphysical 
models hereafter), as they have been discussed in detail 
in Dimmelmeier et al. (2002a,b), and Ott et al. (2007b) 
as well as Dimmelmeier et al. (2007a, b), respectively. (It 
is worth to emphasize, however, the excellent agreement 
found in the hydrodynamical simulations performed with 
three independent numerical codes.) We pay more atten- 
tion instead to the magnetic field evolution. In all our sim- 
ulations an initial magnetic field strength of Bq = 10 10 G 
is considered. This value is an upper limit for the T3M0 
models, since the expected initial toroidal magnetic field 
is of this order (Heger et al. 2005). However, for the D3M0 
models, this field strength is already at (or above) the up- 
per end of the astrophysically expected values. 



For all models we first present results for identical val- 
ues of Bq , in a way that we can study the different effects 
and compare them properly. Afterwards we present the re- 
sults scaled to lower, astrophysically expected values. We 
anticipate that our results can change if some of the several 
assumptions made in our simulations (axisymmetry and 
passive field approximation) are relaxed. An estimation 
and discussion of these effects can be found in Sect. 6.6. 

6.1. Evolution of the magnetic energy parameter 

The evolution of the energy parameter for the magnetic 
field /3 mag = -Emag/l-E-bl can be seen in Fig. 2 for model 
A1B3G5 of our sample. In order to analyze the amplifica- 
tion of the magnetic field, we separate the effects of the 
different components of the magnetic field into (3 V for the 
toroidal component and /3 po i = /3 m ag — Pip for the poloidal 
component, which are also plotted in the figure. As the 
collapse proceeds the magnetic field grows by at least two 
reasons: First, the radial flow compresses the magnetic 
field lines, amplifying the existing poloidal and toroidal 
magnetic field components. Second, during the collapse of 
a rotating star differential rotation is produced and in- 
creased, even for rigidly rotating initial models (see e.g. 
Dimmelmeier et al. (2002b)). Hence, if a seed poloidal field 
exists, the fi-dynamo mechanism winds up the poloidal 
field lines into a toroidal component. This (linear) ampli- 
fication process generates a toroidal magnetic field com- 
ponent, even from purely poloidal initial configurations. 
The toroidal component of the magnetic field is affected 
by the two effects, while the poloidal field is only ampli- 
fied by radial compression of the field lines. Thus, even if 
the initial magnetic field configuration is purely poloidal, 
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Fig. 2. Time evolution of the magnetic energy parameters /3 mag 
(solid line) , f) v (dashed line) , and /3 po i (dashed-dotted line) for 
models A1B3G5-D3M0 (top panel) and A1B3G5-T3M0 (bot- 
tom panel). 

the toroidal component dominates after some dynamical 
time. To study the differences in the evolution of the mag- 
netic field depending on the initial magnetic field we now 
describe in detail the features of model A1B3G5 with dif- 
ferent initial magnetic field configurations. 

In model A1B3G5-D3M0 the initial magnetic field is 
entirely poloidal. The top panel of Fig. 2 shows that V 
(dashed line) grows much faster than /3 po io, particularly 
after bounce (t ~ 30 ms) when the radial compression 
mechanism stops. We note that the magnetic field consid- 
ered is weak enough not to affect the dynamics, with the 
final /3 mag < 1. 

If we consider a purely toroidal magnetic field initially, 
as model A1B3G5-T3M0, the only amplification mecha- 
nism present in our simulations is the radial compression, 
since no poloidal field can be wound up. The bottom panel 
of Fig. 2 shows the behaviour of /3 ma g for model A1B3G3- 
T3M0. It is important to notice that during the collapse 
Pma g hardly grows (for other models of the T3M0 series 
it even decreases) since the radial compression is a very 
inefficient mechanism to amplify the magnetic field. As 
a result, for some models the final PNS is "less magne- 
tized" than the progenitor core in the sense that /3 mag at 
bounce is smaller than it is before the collapse. We note 
that the evolution of this kind of purely toroidal mod- 
els could change completely if the axisymmetry condition 
were removed, since in three dimensions there are mecha- 
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Fig. 3. Time evolution of the magnetic energy parameter /3 mag . 
The top panel shows /3 mag for model A1B3G5-DT3M0 (solid 
line) and the comparison with the composition (dot-dashed 
line) of models A1B3G5-D3M0 (dashed line) and A1B3G5- 
T3M0 (dotted line). The bottom panel shows the evolution 
of /3ma g for the same hydrodynamic model (A1B3G5) with an 
initial value of the toroidal field of = 10 10 G, and varying 
values of the initial poloidal field S; o ° lo , from 10 6 to 10 10 G. 



nisms that can transform a toroidal magnetic field into a 
poloidal one. Some of these mechanisms are discussed in 
Sect. 6.6 below. 

To check whether the "composition rule" (see 
Sect. 2.3) is valid we consider next a mixed configuration 
of poloidal and toroidal magnetic fields at the beginning of 
the simulation (model A1B3G5-DT3M0). The top panel 
of Fig. 3 shows with a solid line the time evolution of 
/3 mag for model A1B3G5-DT3M0 and with a dot-dashed 
line the composition of the individual values for /3 mag in 
models A1B3G5-D3M0 and A1B3G5-T3M0 with identi- 
cal initial field strengths. (The separate evolutions for the 
latter are also included in the plot as dashed and dotted 
lines, respectively.) The agreement of the two evolution 
paths is remarkable, which shows that the "composition 
rule" works properly for our simulations. Therefore, we 
can use it to obtain any desirable composition of mag- 
netic fields with a single hydrodynamic evolution of the 
two models D3M0 and T3M0. For the particular composi- 
tion showed in this model, the final value of /3 mag depends 
very weakly on the initial toroidal magnetic field compo- 
nent. In other words, the structure of the magnetic field of 
the PNS will depend almost exclusively on the radial com- 



12 



Cerda-Duran et al.: Relativistic passive- magneto-rotational core collapse 




t - t [ms] 

Fig. 4. Time evolution of the magnetic energy parameter /3 ma g 
for all simulated models with initial purely poloidal magnetic 
field (label D3M0). 

pression of the initial poloidal component of the magnetic 
field. 

Next, we consider a "composition" of these models 
with different initial magnetic field strength. We keep 
the initial toroidal component fixed at a realistic value, 
B* = 10 10 G, and change the initial poloidal component 
in a range that spans from B*® lo — 10 10 G down to the 
astrophysically more realistic value of 10 6 G. The bottom 
panel of Fig. 3 shows the time evolution of /3 m ag for these 



different configurations. For lower values of B*^ Q , 



the 



f2-dynamo mechanism becomes increasingly slower and 
the initial toroidal component becomes important for the 
magnetic field configuration of the PNS. In the lowest ini- 
tial poloidal field case analyzed, the magnetic field of the 
PNS is completely toroidal and depends exclusively on the 
initial magnetic field configuration. 

The remaining computed models of our sample, includ- 
ing those with microphysics, behave qualitatively in a very 
similar manner, although quantitative differences can be 
found in the amplification of the magnetic field during the 
collapse, and the amplification rates after bounce due to 
the fi-dynamo. Fig. 4 shows the evolution of the magnetic 
energy parameter /3 mag for all the simulated models with 
initial purely poloidal magnetic field (label D3M0). For 
all models we find the following relation between the col- 
lapse time and the amplification rate of the magnetic field 
after bounce (which, however, does not hold for model 
A2B4G1-D3M0, as this is the only model of our sample 
for which the collapse is halted not by the stiffening of the 
EOS, but rather by centrifugal forces at subnuclear den- 
sities; cf. Dimmelmeier et al. 2002b): Models with large 
collapse times, such as all microphysical models as well 
as the simplified model A1B3G3-D3M0, exhibit a more 
efficient amplification of the magnetic field as compared 
to the rapid collapse models (G5 series). To quantify the 
differences between the models we estimate the time scale 
tq for the amplification of the magnetic field by fitting the 
post-bounce evolution of /3 m ag to 



A, 



The resulting values can be found in Table 2. The time 
scale should depend on the central angular velocity f2 c of 
the PNS, and on the strength of the poloidal magnetic 
field that can be wound up, which can be estimated from 
/Spoio. Hence, the following expression should be valid in 
the most efficient scenario (see Appendix B for details): 



Til = 



polo 



(35) 



To check this relation we plot in the top panel of Fig. 5 the 
value of the fit for to versus the value from the previous 
analytic expression. Apparently for all models the growth 
time of the Sl-dynamo is always larger than that of the 
most efficient situation (solid line in the figure), and cor- 
responds to a fraction (30% -90%) of the upper limit (35). 
This relation shows that in order to obtain higher amplifi- 
cation rates of the magnetic field not only strong rotation 
is needed, but also a sufficient compression of the poloidal 
magnetic field during the collapse. 

Furthermore, we also find a relation between the value 
of /3 P oio and the mass enclosed in the neutrino sphere 3 , 
Mpns hereafter (see bottom panel of Fig. 5). Since most 
of the magnetic field lines compressed by the collapse are 
located inside the neutrino sphere, it is easy to understand 
that more massive PNS have higher magnetic energies. 
The fit to a power law of the data shown in Fig. 5 yields 



/3 polo = (3.2 ± 0.5) x 

A/pNS 



10" 



Mr. 



(1.6±0.2) 



B* 



10 10 G 



(36) 



As discussed in detail by Dimmelmeier et al. (2007a, b), in 
the microphysical models the mass of the homologously 
collapsing inner core at bounce has a value of ~ 0.5 Mq 
(for the rotation rates considered here). This is also con- 
sistent with the high mass -Mpns ~ 0.8 Mq of the PNS in 
these models, as shown in Fig. 5. To obtain masses in this 
range in models with a simple matter treatment, the adia- 
batic index would require a value 71 > 1.32, which is close 
to 4/3. Already for moderate rotation, this choice would 
cause the core to undergo multiple centrifugal bounces at 
densities lower than nuclear density (as exemplified here in 
model A2B4G1), which is a dynamical behavior that does 
not occur at all in microphysical models (Dimmelmeier 
et al. 2007a,b, see also the related discussion in Sect. 6.5). 
Therefore, only the microphysical models feature a col- 
lapse to a PNS that has both high densities and is in ad- 
dition comparably heavy. This combination, which cannot 
be realized with the simplified models, explains the higher 
growth rates of the magnetic field due to the 17-dynamo 
observed if improved microphysics is taken into account. 

Combining Eqs. (35) and (36) we can establish an up- 
per limit to the growth rate of the magnetic field due to the 



(34) 



3 In all models, we define the neutrino sphere as the sur- 
face inside the core where the density equals the trapping den- 
sity ptrap ~ 2 x 10 12 g cm -3 . In the microphysical models, 
above this density neutrinos are assumed to be trapped in the 
medium (Liebendorfer 2005). 



Cerda-Duran et al.: Relativistic passive-magneto-rotational core collapse 



13 



100 - 



li 3 10: 



10 



100 



1/2 



2Q. B , [s] 

c 'polo L J 




M pNS /M 

Fig. 5. Top panel: relation between the fitted values for the 
growth time to of the fi-dynamo and the upper limit given 
by 2 f^ < r 1 /3 po io 2 • The solid line represents the upper limit, the 
dashed line corresponds to 50% of that limit, and the dotted 
line to 10%. Bottom panel: magnetic energy parameter f3 po \ of 
the poloidal component of the magnetic field at bounce versus 
the average mass Mpns of the PNS after bounce, and its fit 
to a power law. Error bars are shown for models with mass 
variations larger than 5% after bounce. Both microphysical 
models (filled circles) and simplified models (open circles) are 
presented. Newtonian models (crosses) are also shown, but are 
not used in the fit. 



17-dynamo using only hydrodynamic quantities, namely 
fi c and Mpns, an d the strength of the magnetic field in 
the progenitor, Bq. This limit is given by 



(11.18 ±0.9) x 
1 ms _1 \ ( M G 



fir 



Mi 



PNS 



(0.8±0.1) 



10 10 G 



s. 



(37) 



This relation can be very useful to estimate how fast the 
magnetic field grows in a collapsed star, under the assump- 
tion of a weak magnetic field and with a similar poloidal 
configuration in the progenitor, using data from purely 
hydrodynamical simulations (with no magnetic fields). As 
a proof of consistency and in order to assess the quality of 
this estimate we have computed tq with this method. We 
find that in all cases the estimate is a lower limit for the 
actual value of rn obtained from the numerical simulations 
and deviates by at most 30%. 



6.2. Convection 

One of the most important features that can affect the 
evolution of the magnetic field in stellar core collapse to 
a PNS is the presence of convection. We present here a 
detailed analysis of this effect in our simulations. Since 
in all of our models the magnetic field is weak, the dis- 
cussion can be performed without considering its influ- 
ence. We also note that due to the approximations made 
in our simulations, specifically the lack of a consistent neu- 
trino transport scheme, our findings regarding convection 
should not be considered as definite. 

The stability conditions for a rotating star are given 
by the so-called Solberg-H0iland criteria (Tassoul 1978), 



Cshi =g-B + J- Vw > 0, 
C SH2 = (g x VtJ7)(B x J) > 0, 



(38) 



where g is the gravitational acceleration, and the buoy- 
ancy and rotational terms are respectively given by 



B = 



Vp VP 
P PT l 



J = i-V(fi 2 tu 4 ), 



(39) 



with Ti = (9 In P '/ '9 In p) S; y e . Note that in the first con- 
dition of Eq. (38), iV 2 = g ■ B is the Brunt-Vaisala, fre- 
quency and k 2 — J ■ Vtu is the epicyclic frequency. In the 
case of either no rotation or uniform rotation the Solberg- 
H0iland criteria reduce to the well known Schwarzschild 
criterion, N 2 > 0. If one of the two conditions is not sat- 
isfied, convective instability develops. Following Miralles 
et al. (2004), the time scale of the fastest growing mode 
can be computed as 



t S h = Im 



C S 



HI 



2\/ C shi 



4C: 



SH2 



-1/2' 



(40) 



It is very useful to express the buoyancy terms in the 
conditions (38) in terms of the contributions of the entropy 
and electron fraction gradients, 



B = £Vs + SVY e . : 



(41) 



with £ = -d\nP/ds\ p , Y jTi and 5 = -dlnP/dY^/Y^ 
We point out that the Solberg-H0iland criteria are valid 
exactly only in Newtonian gravity, and thus we use them 
here only as estimates. In order to assess the influence of 
general relativistic corrections, we also evaluate Eq. (38) 
using covariant derivatives with respect to the CFC met- 
ric, which yields very similar results. Note also that the 
Solberg-H0iland criteria are based on a local instability 
analysis, while the convection observed in our simulations 
covers extended regions. 

In Fig. 6 we show the extent of the convectively unsta- 
ble regions according to the Solberg-H0iland criteria (38) 
after core bounce for models of the series s20AlB5, by 
plotting the time evolution of angle-averaged values for 
the convective growth time scale tsh- From this figure it 
becomes apparent that two regions are susceptible to de- 
veloping instabilities: the region just below the neutrino 
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Fig. 6. Time evolution of the angular averaged value of the 
growth time scale tsh, which indicates convectively unstable 
regions according to the Solberg-H0iland criteria. White re- 
gions represent convectively stable regions. The top half of the 
figure shows angular averages near the pole (0 < 6 < 7r/6), 
while the bottom half shows the corresponding averages near 
the equator (n/3 < 6 < 7r/2). The shock radius (dashed lines), 
neutrino sphere radius (dotted line), and radius of shock for- 
mation (dash-dotted line) are also displayed. 



sphere (between about 20 km and 40 km) and extended 
regions behind the shock. The innermost 2 km of the star 
are also convectively unstable, but we suspect that the 
small negative entropy gradient responsible of this unsta- 
ble region is a numerical artifact of the inner boundary, re- 
lated to the so-called wall heating effect commonly appear- 
ing in shock reflection experiments (Donat & Marquina 
1996). In our simulations of models of the series S20A1B5, 
convective motions indeed occur in those unstable regions 
as predicted by the instability criteria, as well as in the 
surrounding regions due to overshooting. We also find that 
the time scale of the onset of the observed instability is 
correctly estimated by Eq. (40). 

Below the neutrino sphere (20-40 km), convection sets 
in inmediately after bounce, with typical maximum veloc- 



ities of about 2 x 10 4 km s _1 . The velocities progressively 
decrease until the end of the simulation (at about 65 ms 
after bounce) with average values around 100 km s~ , al- 
though convection does not disappear completely. Behind 
the shock (100-200 km), the typical convective velocities 
are of the order of 1000 km s - , with maximum values in 
some regions of 10 4 km s~ . This magnitude remains until 
the end of the simulation. 

For a more detailed analysis we separately evaluate 
the different contributions in the Solberg-H0iland crite- 
ria (38) with B in the form of Eq. (41). Since the radial 
gradient of Y e is positive (as deleptonization is stronger 
towards the center during the collapse), this has an sta- 
bilizing effect against convection. Similarly, rotation also 
lM, suppresses convection, since the epicyclic frequency n 2 is 
^ positive everywhere. Convective instability can thus only 
appear in regions with a sufficiently large negative radial 
entropy gradient. Such a gradient occurs in the region al- 
ready swept by the shock front. Shock heating creates en- 
tropy most strongly close to the neutrino sphere at a ra- 
dius of about 30 km (see Fig. 6), producing a steep gradi- 
ent there 1 ms after core bounce. Behind the shock front, 
which then propagates to larger radii at lower densities 
and decelerates, another region with a negative gradient 
also appears. All our microphysical models show very sim- 
ilar qualitative behavior with some variations due to dif- 
ferent angular momentum distribution and the description 
of matter. 

In models with slower rotation (i.e. the S20A1B1 se- 
ries), strong convection sets in immediately after the oc- 
curence of the negative entropy gradient close to the neu- 
trino sphere. For models with very little rotation (which 
have not been considered in this work), such convective 
overturn is strong enough to be clearly visible in the 
post-bounce gravitational wave signal (Dimmclmcicr et al. 
2007a, b). Within about 20 ms after core bounce, con- 
vection has managed to smooth out the entropy gradi- 
ent around the neutrino sphere, thus removing the con- 
dition for sustained convection. Accordingly, convection 
is strongly damped, the vortices disappear quickly, and 
the low-frequency contribution to the gravitational wave 
signal is no longer visible. This fast convective transient 
near the neutrino sphere has been observed in numeri- 
cal simulations without any neutrino treatment (see e.g. 
Burrows & Fryxell 1992; Miiller & Janka 1997), and also 
in simulations using a neutrino diffusion scheme (Swesty 
& Myra 2005), although in the latter case the time scale 
for damping of convection is shorter (~ 10 ms) than in 
our case. However, in simulations including state-of-the- 
art Boltzmann neutrino transport (Miiller et al. 2004), a 
few ms after core bounce no significant convection remains 
in this region, and no traces in the gravitational wave sig- 
nal can be found. We attribute this disagreement with our 
results to the simplified neutrino treatment in our models, 
which cannot properly take into account the deleptoniza- 
tion of the PNS after core bounce. As the deleptonization 
of the PNS is initially strongest when the shock travels 
through the neutrino sphere, we expect the most signifi- 
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cant inaccuracies of our formulation there. We therefore 
conclude that the convection over <~ 20 ms, which we ob- 
serve in the neutrino sphere region, is an artifact that 
should disappear once a more realistic neutrino descrip- 
tion is included. 

In more rapidly rotating models, the stabilizing ef- 
fect of rotation in the Solberg-Efoiland criteria prevents 
the strong transient we find in the slowly rotating mod- 
els from developing, and significantly weaker convection 
is present in this region. However, irrespective of rota- 
tion, convection vortices are formed behind the decelerat- 
ing shock front. On post-bounce evolution times of several 
10 ms, the weak but persistent convection is unable to re- 
move the entropy gradient behind the shock, except near 
the rotation axis, where the specific angular momentum 
is smaller, and convection is stronger. 

Rotation also influences the shape of the convective 
cells. If the buoyancy terms in the Solberg-tfoiland crite- 
ria (38) are much larger than the rotation terms, the con- 
vective cells show no preferred direction. We observe this 
feature particularly in models with slower rotation (the 
S20A1B1 series), and to a lesser degree also in other con- 
vectively unstable models in the first few ms after bounce. 
If the buoyancy terms are comparable in magnitude to the 
rotation terms, convection develops preferredly parallel to 
the rotation axis (see e.g. Miralles et al. 2004). This effect 
is present in our microphysical models at later phases, as 
the entropy gradient has already been partially smoothed 
out and the buoyancy terms have become smaller. 

In contrast to the microphysical models, which show 
remarkable convection in the PNS and behind the shock 
front, models with a simplified matter treatment exhibit 
either no convection at all, or only close to the neutrino 
sphere (in the case of models of the A1B3G5 series). This is 
a consequence of using the hybrid EOS in the latter mod- 
els, which is unable to properly decelerate the shock after 
core bounce and turn it into an accretion shock. Hence 
in these models the entropy gradient is mostly positive 
behind the shock. 

6.3. Structure of the magnetic field 

The main qualitative differences between the various mod- 
els become apparent when we study the detailed struc- 
ture of the magnetic field of the resulting PNS. In Fig. 7 
we show two-dimensional snapshots of selected hydrody- 
namic and magnetic field variables at the final time of 
the simulations for two representative models of our sam- 
ple, namely model A1B3G5-D3M0 (top panels) and model 
s20AlBl-D3M0 (bottom panels). For typical simulations 
with initial poloidal magnetic fields (D3M0 models) the 
resulting PNS has two clearly distinct parts (see left pan- 
els of Fig. 7) : an inner region with a size of ~ 10 km, where 
nuclear density is exceeded and which is almost rigidly ro- 
tating, and a surrounding shell extending to the neutrino 
sphere at <~ 30 km, with subnuclear densities and which 
is strongly differentially rotating. These two parts are also 



visible in the distribution of the magnetic field (see center 
and right panels of Fig. 7). The inner region has a mixed 
toroidal and poloidal magnetic field configuration, with 
both components having similar strength, which results 
in a hclicoidal structure aligned with the rotation axis. As 
this part of the PNS is almost rigidly rotating and practi- 
cally in equilibrium, the magnetic field hardly evolves in 
time. On the other hand, the outer shell is differentially 
rotating; thus the toroidal magnetic field component dom- 
inates and grows linearly with time due to the Sl-dynamo 
mechanism. 

If we compare the microphysical with the simplified 
simulations, we find that some significant morphological 
differences arise due to the stronger convection in the mi- 
crophysical models just below the neutrino sphere. These 
motions affect the magnetic field, since they twist the 
poloidal magnetic field lines, generating a much more com- 
plicated structure of the poloidal field for those models. 
In particular those strong meridional currents distort the 
magnetic field in such a way that in some regions the 
poloidal component changes direction with respect to the 
rotation axis (see e.g. bottom-right panel of Fig. 7). This 
produces a negative effect in the fi-dynamo as in these 
regions the magnetic field is wound up in the opposite di- 
rection. However, the overall f2-dynamo mechanism seems 
not to be affected in a significant way by these local effects. 

Model A4B5G5-D3M0 has to be discussed separately, 
since it has initially significantly stronger differential rota- 
tion and more angular momentum than the other models. 
As a result this model undergoes a core bounce due to 
centrifugal hang-up before reaching nuclear density. Its 
structure is toroidal with an off-center maximum density. 
Although it exhibits stronger differential rotation at the 
beginning compared to the other models, and the amplifi- 
cation process during collapse is thus more efficient, after 
bounce its angular velocity f2 is smaller (as the PNS is 
less compact) and therefore the linear amplification due 
to O-dynamo is less pronounced. The main differences in 
the magnetic field structure of its PNS with respect to 
the other models are that, first, the fi-dynamo is active 
not only in the high-density torus, but also in the cen- 
tral lower-density region, and, second, the strong merid- 
ional currents twist the magnetic field lines around the 
torus. However, we point out that in the investigated 
range of initial rotation configurations all microphysical 
models are significantly less influenced by rotation than 
the simplified models (like A4B5G5), and that even for 
rather extreme rotation such collapse dynamics, leading 
to a toroidal structure, is strongly suppressed if an ad- 
vanced description of microphysics is used (which is in 
accordance with the comprehensive parameter study by 
Dimmclmeier et al. 2007a). 

In the models with initially purely toroidal field at the 
beginning (T3M0 series), a poloidal field cannot emerge in 
axisymmetry. Hence, the final magnetic field structure of 
the PNS consists of a stationary and entirely toroidal mag- 
netic configuration with the highest field strengths found 
in the high density regions. As the rotational profile does 
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Fig. 7. Configuration of the innermost region of the collapsed star at the end of the evolution for models A1B3G5-D3M0 (top 
panels; t — 60 ms) and S20A1B1-D3M0 (bottom panels; t — 142.5 ms). The left panels show the rest mass density as logp in 
units of g cm -3 , overplotted by the meridional velocity field (v r ,v e ) (arrows), and isocontours of the specific internal energy 
e. The center panels display the logarithm of the poloidal component of the magnetic field, log |B po io|, in units of G and the 
magnetic field lines in the r—6 plane. The right panels show B v /\B po i \. All axes are in units of km. 



not affect the distribution of the magnetic field, the dif- 
ferent regions of the PNS are not visible in the structure 
of the magnetic field. 



6.4. Comparison with Newtonian results 

In order to study the general relativistic effects in the evo- 
lution of the magnetic field, we choose a subset of our 
simulations with the hybrid EOS to represent the rela- 
tivistic version of some of the models of Obergaulinger 
et al. (2006a, b). Their first paper is devoted to Newtonian 
simulations of magneto-rotational core collapse, while in 
their second paper an effective relativistic gravitational 
potential was used to mimic general relativistic effects 
(while still keeping a Newtonian framework for the hy- 
drodynamics; TOV models in their notation). Since in 
contrast to their work we use the passive field approxi- 
mation, the comparison can only be made with the low 
magnetic field models presented in that work, namely the 
"M10" models. In these models the magnetic field does 
not affect the collapse dynamics and our approximation 
is valid. Although there are no qualitative differences be- 
tween Newtonian and general relativistic models (aside 
from those coming purely from the hydrodynamics as de- 
scribed in Dimmclmeier ct al. 2002a, b), some dissimilari- 
ties can be found in the magnetic field strength and am- 
plification rates after core bounce. 



We have studied the evolution of the magnetic energy 
parameter /3 ma g for the various models, and plot the re- 
sults in Fig. 8. Note that for the same initial magnetic 
field, the magnetic field contribution to the magnetic en- 
ergy parameter differs between a purely Newtonian treat- 
ment, a Newtonian formulation with the effective rela- 
tivistic TOV potential, and general relativity. As a conse- 
quence, the initial value of /3 m ag is not the same in these 
three cases. In order to be able to make an unambiguous 
comparison, we scale the magnetic fields such that /3 ma g in 
the initial model is equal to the value in general relativity. 
In general, for a similar hydrodynamic behavior (models 
A1B3G3-D3M0 and A3B3G5-D3M0) the magnetic energy 
attained during the evolution is smaller in the general rel- 
ativistic case than in the Newtonian case (with either reg- 
ular or effective relativistic TOV potential). 

The winding up of magnetic field lines is the main 
mechanism responsible for the increase of the magnetic 
field during the collapse. Therefore the amplification rate 
for /3 mag is determined by what rotation rate is reached 
and also by how strongly the poloidal component of the 
magnetic field is compressed. In the general relativistic 
case both higher densities and also stronger rotation are 
achieved (Dimmelmeier et al. 2002b). To investigate the 
impact of general relativistic gravity on the magnetic field 
compression, we consider /3 po io as this quantity is the seed 
for the fi-dynamo. In general relativity the PNS has in 
general a smaller mass -Mpns than in the corresponding 
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Fig. 8. Comparison of the time evolution of the magnetic energy -B mag for models A1B3G3-D3M0 (left panel), A3B3G5-D3M0 
(center panel) and A2B4G1-D3M0 (right panel). The line styles represent simulations performed in general relativity (solid 
lines), a purely Newtonian treatment (dotted lines), and Newtonian hydrodynamics with an effective relativistic TOV potential 
(dashed lines). 



Newtonian simulation of the same model. Following the 
relation established in Sect. 6.1 (see the bottom panel of 
Fig. 5), the smaller PNS mass in the general relativistic 
simulation leads to a lower value of (3 V0 \ . Therefore in that 
case, despite the larger f2 c the much smaller magnitude of 
/3 P oio results in a longer time scale for the f2-dynamo via 
Eq. (35), and hence a smaller growth rate of the magnetic 
field. 

In the multiple centrifugal bounce model A2B4G1- 
D3M0, general relativistic effects lead to a bounce at sig- 
nificantly higher maximum densities than in Newtonian 
gravity. Therefore, this is the only investigated model 
where Mpns, and consequently /3 po io as well as /3 m ag are 
larger in the general relativistic simulation. 

6.5. Gravitational waves 

We calculate the gravitational wave output from all of 
our simulations using the Newtonian quadrupole formula 
given in Eq. (22), which includes the magnetic terms. 
Thus, the quadrupole wave amplitude A^q, which is re- 
lated to the dimensionless quadrupolar strain amplitude 
/ii uad as 



^quad 



l 



15 sin 2 0^1. 

7T R 



(42) 



contains the contribution A^Q m&g corresponding to the 
magnetic field. Here h quad is the only independent com- 
ponent of the radiative part hf^ ad of the spatial metric 
as given by Eq. (20). In order to understand how the 
magnetic field affects the waveforms, we also separately 
compute ^fomag- The resulting waveforms for some rep- 
resentative models are shown in Fig. 9. As the magnetic 
field is very low at all times, b 2 <C p, the component of 
the gravitational wave due to the magnetic field is several 
orders of magnitude smaller than A®q ■ 

Therefore, during the core bounce and the im- 
mediate post-bounce phase, the waveforms we ob- 
tain are practically identical to the ones presented 
for the same model setup in Dimmelmeier et al. 



(2002a) (for the simplified models) and Dimmelmeier 
et al. (2007a) (for the microphysical models), which 
can also be downloaded from a freely accessible wave- 
form catalog at www.mpa-garching.mpg.de/relJrydro/ 
wave_catalog. shtml. The values for A^q lie in the range 
between about 30 cm and 3000 cm, which translates to a 
/i quad of roughly 3 x 10~ 22 to 3 x 10~ 20 (assuming a dis- 
tance R — 10 kpc to the source and optimal orientation 
between the source and the detector). 

We also emphasize that all investigated microphysical 
models yield gravitational wave signals known as Type I 
in the literature, i.e. the waveform exhibits a positive pre- 
bounce rise and then a large negative peak, followed by a 
ring-down. This is to be expected, as recent studies using 
the same hydrodynamical model setup (Ott et al. 2007a; 
Dimmelmeier et al. 2007a) have shown that the inclusion 
of microphysics in stellar core collapse simulations sup- 
presses the other signal types, which were associated to 
multiple centrifugal bounce (Type II signals) or rapid col- 
lapse with a very small mass of the inner core (Type III 
signals). 

After bounce, the star reaches a quasi-equilibrium 
state, and thus, the hydrodynamic component of the wave- 
form decreases. At the same time, for models D3M0, the 
magnetic field grows linearly with time. Such a behav- 
ior in the magnetic field produces an increasing gravita- 
tional wave signal, which grows quadratically with time 
due to the dependence on the magnetic field in Eq. (22). 
However, at the end of the simulation, the magnetic field 
component of the waveform is still negligible in compar- 
ison with the hydrodynamic component. It is expected 
that at later times, as the amplification of the magnetic 
field reaches saturation, the influence of the magnetic field 
on the waveform becomes significant, both due to its ef- 
fect on the dynamics and also due to the contribution of 
the magnetic field to the gravitational radiation itself. We 
note, however, that the effect of the MRI could addition- 
ally lead to noticeable changes in the waveforms, provided 
it were able to efficiently amplify the magnetic field (see 
discussion in Sect. 6.6.2). 
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Fig. 9. Absolute value of the gravitational wave amplitude A?o 
(solid line) for models A1B3G5-D3M0/T3M0/DT3M0 (top 
panel) and S20A1B5-D3M0/T3M0 (bottom panel). For the low 
magnetic field strengths considered, the contribution of the 
magnetic field to the waveform is negligible, and the signals 
of series D3M0, T3M0, and DT3M0 are practically identical 
to the purely hydrodynamic waveform. For clarity, the compo- 
nent Afo mag from the magnetic field is also plotted for mod- 
els A1B3G5-D3M0 (top panel, dashed line), A1B3G5-T3M0 
(top panel, dotted line), A1B3G5-DT3M0 (top panel, dashed- 
dotted line), S20A1B5-D3M0 (bottom panel, dashed line), and 
S20A1B5-T3M0 (bottom panel, dotted line). 



For models T3M0, on the other hand, the component 
of the waveform due to the magnetic field is even smaller 
than for the D3M0 models. This is a consequence of the 
inefficient amplification of the magnetic field via the radial 
compression. After bounce, the magnetic component of 
the waveforms in models T3M0 does not grow, and hence 
it is not expected to dominate the waveform later in the 
evolution, unless other processes amplifying the magnetic 
field were present. 

6.6. Amplification of the magnetic field 

Different mechanisms that amplify the magnetic field can 
act during a core collapse or the subsequent evolution of 
the newly formed PNS. This issue is of great importance, 
since the evolution of the PNS during its first minute of 
life until a cold NS forms can change drastically depend- 
ing on the initial conditions at formation. One of the most 
important aspects is the distribution of angular momen- 



tum. A highly differentially rotating PNS can be subject 
to various types of instabilities, such as the dynamical low- 
(3 instability, the classical bar-mode instability (which is 
unlikely to occur in a PNS on dynamical time scales as it 
requires very high values of /3), or the secular CFS insta- 
bility. Such instabilities are potential sources of detectable 
gravitational radiation. Therefore, a natural question that 
arises is whether the magnetic field is going to grow suffi- 
ciently fast to act on the PNS dynamics by flattening the 
rotation profiles (and therefore preventing the instabili- 
ties to develop), or whether, instead, the growth process 
may take a few seconds, allowing the instabilities to grow 
and the accompanying gravitational waves to become de- 
tectable. A number of effects can amplify the magnetic 
field shortly after PNS formation. In the following, we 
discuss these effects and estimate their importance for our 
models of core collapse 4 . 

6.6.1. fi-dynamo 

Within our passive field approximation we can only com- 
pute the amplification rates for the O-dynamo, for which 
the magnetic field grows linearly with time; therefore /3 m ag 
grows quadratically with time (see Appendix B). The time 
scale tq of this amplification process and the estimated 
time i sat at which the field saturation begins are given in 
Table 2. In the fastest case of our model sample, which 
occurs for model A1B3G3-D3M0, saturation is reached at 
about 300 ms, and in most other cases, the f2-dynamo sat- 
urates at times larger than 0.5 s. Note, however, that these 
estimates depend on the initial magnetic field strength, 
which is chosen to be Bq = 10 10 G. For lower values of 
the magnetic field these time scales can be scaled as (see 
Eq. B.5) 

/10 10 G\ /10 10 G\ , ao . 

TQ ~ TfilO y— -^—J ) ^sat ~ isatlO ^— J • (43) 

We recall that stellar evolution calculations predict that 
in a progenitor core the poloidal component of the mag- 
netic field can initially have a strength of about 10 6 G 
(Hcger et al. 2005). For such an initial magnetic field the 
saturation time scale becomes several hours. This makes 
the fi-dynamo a very inefficient mechanism to amplify the 
magnetic field during core collapse and bounce, unless the 
progenitors are highly magnetized (B > 10 10 G) for which 
the saturation could be reached within a few dynamical 
time scales. The magnetic field at the saturation is inde- 
pendent of the initial magnetic field strength and of the 
order of - 10 16 G. 



4 For these estimates we utilize the Newtonian limit, since 
most of the work on linear analysis of instabilities has not 
yet been extended to general relativity. Furthermore, for an 
approximate assessment, the restriction to a Newtonian treat- 
ment appears sufficiently accurate. 
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6.6.2. Magneto-rotational instability 



be roughly estimated as 5 



There are other magnetic field amplification processes that 
our simulations cannot account for, but for which it is nev- 
ertheless possible to estimate the growth rates. It has been 
suggested that the magneto-rotational instability could 
amplify the magnetic field from arbitrary weak fields up to 
values where equipartition between the magnetic field en- 
ergy and the rotational kinetic energy is reached (Akiyama 
et al. 2003). However, our analysis shows that in the con- 
text of core collapse such an amplification is still an open 
issue. We proceed next to describe the MRI and the un- 
certainties related to its effect on the amplification of the 
magnetic field in core collapse. 



Linear regime: The MRI is a shear instability that gener- 
ates turbulence and results in an amplification of the mag- 
netic field in a differentially rotating magnetized plasma 
(Balbus & Hawley 1991, 1992), redistributing angular mo- 
mentum in the plasma. Linear analysis shows that if the 
magnetic field strength is very low, as in our case, the sta- 
bility criteria for the MRI in the Newtonian limit (Balbus 
1995) are 



T MRI = 47T 



Cmrii = g-B + TZ-Vm>0, 
C M ri2 = {g x Vw)(B x TZ) > 0, 



(44) 



where TZ = m\7(fl 2 ). Note that these criteria are very simi- 
lar to the Solberg-H0iland criteria (38) for convection, but 
use an angular velocity gradient TZ instead of an angular 
momentum gradient J . Since in the core collapse scenario 
TZ < is satisfied almost everywhere, it is important to 
compute the buoyancy terms given by B to estimate the 
onset of the MRI. For regions with B > (i.e with a nega- 
tive entropy gradient that is strong enough to compensate 
the positive electron fraction gradient term in Eq. (41)), 
the first criterion is not fulfilled. Furthermore, for regions 
with B < (i.e. a positive or sufficiently small negative 
entropy gradient), the second criterion is neither satisfied. 
This means that the presence of a adequately strong neg- 
ative entropy gradient (which also leads to convective in- 
stability) enhances the MRI, although a positive entropy 
gradient does not affect the condition for MRI instabil- 
ity. Note that this peculiarity is caused by the negative 
value of TZ, and does not happen in the Solberg-H0iland 
criteria (38) for convection, as J > in that case. If at 
least one of the criteria (44) is not satisfied and a magnetic 
field is present, then fluid and magnetic field perturbations 
grow exponentially in time. Neglecting buoyancy terms, 
the time scale for the fastest growing unstable mode can 



dim 



(45) 



which is independent of the magnetic field configuration 
and strength. Only those modes with a length scale larger 
than a critical wavelength will grow (Balbus & Hawley 
1991). This length scale can roughly be estimated as 
Amri ~ 2ttca/Q, where ca = y B 2 / p is the Alfven speed. 
For the typical values attained in the nascent PNS, in 
which the dominant magnetic field is toroidal, the critical 
length scale is 



Amri ~ 62 



B 



* o 



10 10 G 



1 ms 



10 g cm 



1/2 



(46) 



Note that we have scaled the length scale with the typical 
magnetic field strength B* of the progenitor, and not 
with that of the PNS itself. For the poloidal component 
and realistic values of the initial magnetic field (B*° <~ 
10 6 G) this length scale is reduced by several orders of 
magnitude (Amri P o1o ~ 0.6 cm). In any case, resolving the 
scales needed to simulate the MRI is a challenging problem 
as, in the case of weak magnetic fields, the wavelength of 
the fastest growing mode (which is close to the critical 
length scale) is typically much smaller than the available 
grid resolution. 

Non-linear regime: Linear analysis provides tools to deter- 
mine the onset of the instability and the typical time and 
length scales. However, once the perturbations of the mag- 
netic field reach values comparable to the magnetic field 
itself, linear analysis is no longer valid (although in the 
weak field case the perturbations of the fluid variables are 
still small) . The amplification of the magnetic field due to 
the MRI is therefore a nonlinear effect, and can only be 
studied by means of numerical simulations. The appropri- 
ate numerical approach, due to the smallness of the length 
scales necessary to be resolved, are local simulations of the 
MRI in a shearing box. Numerical simulations of this kind 
in three dimensions have been performed by Hawley et al. 
(1995) in the context of accretion discs. They show that 
if the instability condition of linear analysis is fulfilled, 
then the amplification of the magnetic field proceeds by 
the formation of an axisymmetric channel flow. This is 
well understood, since the linear MRI solution is also a 
solution of the nonlinear axisymmetric MHD equations 
(Goodman & Xu 1994). In the ideal MHD limit, the am- 
plification saturates when nonaxisymmetric perturbations 
destroy the channel flow. It is important to emphasize the 



5 We note that Balbus & Hawley (1991) derived a compli- 
cated expression including bouyancy terms which, however, is 
only valid in the equatorial plane. To the best of our knowl- 
edge the timescale for the fastest growing mode in the general 
case has not been computed so far. It would require solving 
the dispersion relation, a task out of the scope of this paper. 
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necessity of performing three-dimensional simulations in 
the shearing box since, in axisymmetry, the channel flow 
is not destroyed and any magnetic field is able to grow con- 
tinuously, reaching saturation only when the MRI length 
scale is of the order of the region in which the MRI is 
present (Hawley & Balbus 1992). 

For a magnetic field distribution with zero mean at 
large scales, the amplification proceeds from arbitrarily 
weak fields (Hawley et al. 1996) and saturates irrespec- 
tive of the initial magnetic field at average values of 
Pmag/Po ~ 0.01, where Pq is the initial gas pressure. If a 
mean magnetic field is present (as in our case), the satu- 
rated magnetic field depends on the initial magnetic field 
strength. In the most favorable case of a vertical mag- 
netic field, the total amplification by the MRI is only 
about a factor 20 of the original field, this amplifica- 
tion being even smaller in the case of a purely toroidal 
field (Hawley et al. 1995). On the other hand, Sano et al. 
(2004) have suggested that for sufficiently weak magnetic 
fields the saturation level could be independent of the ini- 
tial field and equal to that in the zero-mean case. If this 
were confirmed it would mean that, for the weak magnetic 
field strength present in our magneto-rotational core col- 
lapse models (-P m ag/-Po ~ 10~ 8 in the PNS for progeni- 
tors with B*° = 10 10 G), a magnetic field of - 10 16 G 
could be reached on time scales of tmri- Such a strong 
magnetic field would have a significant effect on the dy- 
namics, similar to that observed in numerical simulations 
with highly magnetized progenitors (Obcrgaulinger et al. 
2006b, a; Shibata et al. 2006). In the opposite case, the 
MRI would fail to considerably amplify the magnetic field, 
and for a purely toroidal field the magnetic field should 
grow only by a factor of about 3 according to Hawley et al. 
(1995). 

The inclusion of more complex physics relevant for the 
core collapse scenario (like radiation, diffusion, or resis- 
tivity) can significantly change the amplification process, 
since in the nonideal case the reconnection of magnetic 
field lines seems to be the dominant effect in the satu- 
ration process of the MRI. In general, these effects act 
towards lowering the values of the saturation; for reasons 
of simplicity we do not consider them in this discussion 
(see Hawley 2005, and references therein for a detailed 
review on this topic). 

Furthermore, it has to be noted that all local simu- 
lations of the MRI have been performed in the context 
of Keplerian accretion discs, and, hence, some of the un- 
derlying physical conditions are not valid in the case of 
a PNS. For example, the typical sound speed c s in those 
simulations is of the order of 10~ 3 . Only the parametric 
study performed by Sano et al. (2004) covers a wider range 
of values of c s ~ 10~ 8 -10~ 2 in units of c. However, the 
sound speed in a PNS is higher, c s ~ 10 _1 . Rotational ve- 
locities and profiles are also very different in a disc and a 
PNS. Therefore, appropriate local simulations of the PNS 
scenario should eventually be performed in order to con- 
firm the growth of the MRI for a weakly magnetized PNS, 



and to infer the magnetic field at which the instability sat- 
urates. 

Our results: As the MRI involves a backreaction of the 
magnetic field onto the dynamics, we cannot study this 
effect in our simulations, as we assume the passive field 
approximation. Furthermore, even with "active" magnetic 
fields, both the resolution needed to resolve the MRI 
length scale (~ 10 m) and the requirement for three- 
dimensional simulations are not affordable with present 
computers. Therefore, we are limited to analyzing whether 
our magnetized collapse models are susceptible to devel- 
oping such an instability according to linear analysis esti- 
mates, leaving aside the issue of saturation, whose uncer- 
tainties need a deeper analysis which is beyond the scope 
of this work. In order to estimate how MRI could change 
our results if it were taken into account properly, we deter- 
mine the regions where the MRI instability criteria (44) 
are not satisfied. Inside these regions we calculate the time 
scale for the fastest growing mode. In Fig. 10 we show 
the results for the models A1B3G5-D3M0 (left panel) and 
S20A1B5-D3M0 (right panel). We note that since the onset 
of the MRI is independent of the magnetic field strength, 
provided that a poloidal component exists, any composi- 
tion of D3M0 and T3M0 models has the same instability 
properties as the D3M0 models. 

Our analysis of all computed models shows that during 
the infall phase the MRI is either not possible or the typi- 
cal time scales involved are much larger (i.e. > 10 s) than 
the duration of the collapse itself. Therefore the instabil- 
ity can affect neither the dynamics nor the magnetic field 
strength in that phase. Around the time of core bounce, 
the angular velocity gradient is larger and the MRI time 
scale becomes dynamical. Almost the entire region be- 
tween the shock formation radius (at <~ 10 km) and the 
shock itself is MRI unstable with time scales of the order 
of <~ 1-10 ms. Note that the innermost part of the PNS 
rotates rigidly, and therefore the MRI unstable region that 
appears in the inner 2 km is possibly a numerical artifact 
caused by the probably unphysical negative entropy gradi- 
ent mentioned in Sect. 6.2. Some differences appear when 
comparing microphysical and simplified models. 

A general feature of the microphysical models is the 
post-bounce appearance of a negative entropy gradient 
(regions with N 2 < 0, see Sect. 6.2). This property is 
much less prominent in the simplified models (except in 
model A1B3G5). Thus, the presence of a such gradient 
in the microphysical models enhances the occurence of 
the MRI behind the shock as compared to the simplified 
models (see Fig. 10), since in these regions the cause for 
the instability is mainly the presence of a negative entropy 
gradient. Around the neutrino sphere the presence or ab- 
sence of a negative entropy gradient does not affect the 
onset of the instability since it is caused by the strong 
negative angular velocity gradient. Therefore, only small 
differences can be found in the latter region between the 
simplified and the microphysical models. 
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Fig. 10. Time evolution of the angular averaged value of the growth time scale tmri of the fastest growing mode of the MRI 
for models A1B3G5-D3M0 (left panel) and S20A1B5 (right panel). White regions are stable to the MRI. The top half panels 
show angular averages of tmri near the pole (0 < 8 < 7r/6), while the bottom half panels show these averages near the 
equator (n/3 < 8 < 7r/2). The shock radius (dashed lines), neutrino sphere radius (dotted line), and radius of shock formation 
(dash-dotted line) are also displayed. 



As a result of this analysis, for collapse progenitors 
with a magnetic field smaller than 10 10 G (hence includ- 
ing astrophysically more relevant initial values of 10 6 G), 
we infer that perturbations of the magnetic field are going 
to grow exponentially on dynamical time scales and will 
reach saturation in the unstable regions mentioned above. 
However, the value of the magnetic field at which satura- 
tion appears is still unknown, which is a key issue in order 
to establish the effects of the MRI, if any, on the dynamics. 
Nevertheless, even if the MRI were unable to considerably 
amplify the magnetic field, it could still play a major role 
at late times during the evolution of the PNS, provided 
other amplification mechanisms were capable to increase 
the magnetic field to significant larger values (see below). 
In such a situation the MRI could have an impact on the 
dynamics by transporting angular momentum outwards 
and driving the PNS towards rigid rotation. 



6.6.3. Dynamo mechanisms 

The wind-up process of the magnetic field (fi-dynamo) 
discussed before is a mechanism that works by transform- 
ing the poloidal magnetic field into a toroidal field and ex- 
tracting energy from differential rotation. In axisymmetry 
this process amplifies the magnetic field linearly with time 
as long as differential rotation exists. If the axisymmetry 
condition is relaxed, however, a number of instabilities of 
the toroidal field can transform the toroidal magnetic field 
back into a poloidal magnetic field. This feedback then 
"closes" the dynamo process. 

The first group of instabilities are those related to con- 
vective unstable regions, neutron-finger instabilities (due 
to a negative composition gradient) and, in general, tur- 
bulence. In these cases the cv-effect is the one which closes 
the dynamo in the a-f2-dynamo (Thompson & Duncan 
1993). Computations of this effect (Bonanno et al. 2005) 
suggest that even for a rapidly rotating PNS with a pe- 
riod around 1 ms (i.e. comparable to the models presented 
here), the time scale for the growth of the magnetic field 
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is <~ 1 s. Therefore, this effect is probably not important 
after core bounce on dynamical time scales. However, for 
larger time scales (i.e. several seconds), and if the MRI is 
not efficient enough, this mechanism will most likely am- 
plify the magnetic field, leading to magnetic braking of 
the PNS within a few seconds. 

There are also types of instabilities that can act in 
stably stratified regions, i.e. regions which are convectively 
stable. Spruit (1999) has proposed the Tayler instability 
(Tayler 1973) as a mechanism to close the dynamo. This 
dynamo has been confirmed in numerical simulations by 
Braithwaite (2006a, b). The condition for this kink- type 
instability to grow in the rotating case (to = 1 mode) is 
(Spruit 1999) 

dglnBl sin (9 cos (9 > 0, (47) 

which is satisfied almost everywhere inside the star in our 
simulations. The growth rate of the instability is of the 
order of the Alfven time scale, 

r T = ^ (n«n A ), TT = -^~ ( fi »^A), (48) 

where J1a = cp^/R and R is the typical size of the re- 
gion considered. In case this instability appears, it de- 
stroys the toroidal magnetic field by transforming it into 
a poloidal field which feeds back the amplification of the 
toroidal magnetic field via the fi-dynamo. Therefore, the 
dynamo is only effective if the O-dynamo is able to gener- 
ate a toroidal magnetic field faster than the Tayler insta- 
bility destroys that field, i.e. tt 3> tq. Saturation is then 
reached as tt ~ tq. Note that depending on the system, 
the saturated magnetic field can be weak enough not to 
affect the dynamics. 

If we consider the typical toroidal magnetic field at 
bounce to be 10 13 G (as in the T3M0 models) with a typ- 
ical density in the PNS of p ~ 2 x 10 14 g cm~ 3 and a typ- 
ical size of the inner region of R <~ 10 km, then the time 
scale for the growth of the Tayler instability is strongly 
increased by rotation, 

^ 3 hH l^s-rj hr < (49) 

which we obtain from the il 3> Oa limit of Eq. (48). 

This means that for a typical progenitor with a toroidal 
magnetic field of £?*° ~ 10 10 G, the instability proposed 
by Spruit (1999) is going to be very inefficient in ampli- 
fying the magnetic field. However, on a longer time scale, 
when other mechanisms could amplify the magnetic field 
(e.g. the a-fi-dynamo) , the Tayler instability could also 
become important. 

7. Conclusions 

In this paper we have presented numerical simulations of 
the collapse of rotating magnetized stellar cores in the 
CFC approximation of general relativity, as well as tests 



assessing our numerical approach for solving the ideal gen- 
eral relativistic magneto-hydrodynamics (GRMHD) equa- 
tions. 

As initial models we have set up (either fully or nearly) 
stationary configurations of weakly magnetized stars in 
general relativity, with either toroidal or poloidal (or both) 
magnetic field components. We have used the "test" pas- 
sive field approximation for evolving these initial models, 
for which the magnetic pressure in all cases considered is 
several orders of magnitude smaller than the fluid pres- 
sure. 

We have performed tests to check the accuracy and 
convergence properties for the GRMHD extension of our 
code. For magnetic field quantities we have found an order 
of convergence above 1 in all of the performed tests. These 
results are consistent with the second-order accuracy (in 
space and time) of our numerical scheme, reduced to first 
order only at shocks and local extrema. The errors in all 
of the cases in which the theoretical solution is known are 
below 0.1%, except at shocks, which are correctly captured 
within only few numerical cells thanks to the use of high- 
resolution shock-capturing schemes. 

For the simulations of magnetized core collapse, we 
have considered cases with magnetic fields which are ini- 
tially either purely poloidal (series D3M0) , purely toroidal 
(series T3M0), or a combination of both. The D3M0 mod- 
els are a general relativistic extension of a subset of the 
cases evolved in fully coupled MHD by Obergaulinger 
et al. (2006a, b), who used a Newtonian formulation (ap- 
proximating general relativistic effects to some extent in 
the latter work). One of our aims has been to compare 
the dynamics and gravitational waveforms with their re- 
sults. No qualitative differences have been found in the 
models studied, while quantitatively the strength of the 
magnetic field at bounce and after the collapse are consis- 
tently smaller in general relativity. 

We have also compared simulations of models with 
improved microphysics (employing a tabulated non-zero 
temperature equation of state (EOS) and an approxi- 
mate but effective deleptonization scheme) with the sim- 
ple (though still widely used) analytic hybrid EOS. The 
results show that the microphysical models (i) lead to 
a more complex structure of the poloidal magnetic field 
due to convective motions surrounding the inner region of 
the PNS, and (ii) exhibit a wind-up of the magnetic field 
(il-dynamo) that is more efficient than in the simplified 
models for comparable rotation rates, which is due to the 
larger compression of the poloidal component during the 
collapse. 

We have found a unified explanation for the mag- 
netic energy of all models, independent of the descrip- 
tion of gravity (general relativity or Newtonian) or the 
EOS, which relates the angular velocity and mass of the 
PNS with its magnetic energy and the growth rate of the 
magnetic field due to the fi-dynamo. This relation shows 
that higher rotation rates and masses of the PNS lead 
to stronger magnetic fields. We have shown that it is not 
possible to mimic the conditions of the microphysical sim- 
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ulations using a simplified EOS. Simplified models with a 
mass of the homologously collapsing inner core during con- 
traction and a mass of the PNS after bounce similar to the 
respective masses of the microphysical models (and identi- 
cal initial rotation profiles) will undergo multiple centrifu- 
gal bounces, a behavior that has recently shown to be an 
artifact of the neglect of microphysics (Dimmelmcicr ct al. 
2007a). 

Further differences appear in the appearance of con- 
vective motion in the PNS and behind the shock. This 
convection is more active in microphysical models than in 
simplified ones. In models with slow rotation, strong con- 
vection in the PNS occurs as a transient and disappears 
within a few ten ms after bounce. Evidently, this tran- 
sient is an artifact as it does not appear in simulations 
of similar models with comparable microphysics but using 
Boltzmann neutrino transport (Miiller et al. 2004) instead 
of our simple advection scheme for the electron fraction 
after core bounce. In rapidly rotating models convection 
is not entirely suppressed by rotation but develops and 
persists on longer time scales, albeit at a lower intensity. 

As we have adopted the passive field approximation 
and the investigated magnetic fields are weak in all phases 
of the collapse, the waveforms of the gravitational radia- 
tion emitted by all our models are practically identical to 
the corresponding ones in a purely hydrodynamical simu- 
lation (Ott et al. 2007a; Dimmelmeier et al. 2007a), with 
the contribution due to magnetic fields being several or- 
ders of magnitude smaller than the total signal ampli- 
tude. However, if the MRI could become dominant for 
the dynamics in the post-bounce phase, in a fully coupled 
GRMHD simulation we would expect a clear imprint of 
such an instability on the signal waveform. As expected, 
for the microphysical models we obtain gravitational wave 
signals exclusively of Type I, as all other waveform types 
(in particular the Type II signals associated with multi- 
ple centrifugal bounces) are suppressed if more realistic 
microphysics is taken into account. 

For an astrophysically expected strength of the mag- 
netic field (Heger et al. 2005), where the initial toroidal 
component is much larger than the poloidal one, we have 
obtained a topology of the magnetic field in the PNS 
that is purely toroidal due to the radial compression of 
the initial toroidal component. In this case the time scale 
for the fi-dynamo is very long (several hours). For pro- 
genitors with stronger poloidal magnetic fields, we have 
found that a core/shell structure is formed. Inside the core, 
where nuclear density is exceeded, a mixed configuration 
of a poloidal and a toroidal magnetic field yields a heli- 
coidal configuration of the field lines. In the surrounding 
shell (which extends several 10 km) the poloidal magnetic 
field lines arc wound up due to differential rotation (fl- 
dynamo) , and shortly after core bounce the magnetic field 
is dominated by the toroidal component. The growth time 
scale for the toroidal component due to this process is, in 
the best case scenario, several 100 ms. 

We have also estimated the growth times for several 
other instabilities that could appear if the passive field 



approximation or the restriction to axisymmetry are re- 
moved. Among these the MRI is apparently the fastest 
growing instability, although it remains unclear if it is go- 
ing to amplify the magnetic field sufficiently (from the 
initially weak field values) to affect the dynamics at all. In 
addition, we have found that the inclusion of microphysics 
could enhance the MRI, since the regions behind the shock 
exhibit a negative entropy gradient, resulting in a growth 
time of ~ 10 ms for the MRI. However, the influence of 
our simplified neutrino treatment or the effects of an al- 
ternative microphysical EOS must still be investigated in 
detail. 

In the event that the MRI were unable to sufficiently 
amplify the magnetic field in the PNS (which is still 
an open issue), the main amplification mechanism would 
probably be the a-O-dynamo, which can amplify the mag- 
netic field to values where the magnetic energy is in 
equipartition with the rotational kinetic energy on a time 
scale of, at least, several seconds. The study of this ef- 
fect is well beyond the goals of the work presented in 
this paper, since the required time scales are much longer 
than those affordable with current numerical magneto- 
hydrodynamical (MHD) codes. Moreover, the underly- 
ing physics necessary to be included (like neutrino trans- 
port, diffusion, radiation, and cooling) is far more com- 
plex. However, in the light of the results presented here, 
in which astrophysically expected values for the magnetic 
field have been adopted, we can speculate about the fol- 
lowing scenario. If the MRI is ineffective, after core bounce 
the magnetic field does not grow significantly strong dur- 
ing one (or maybe several) seconds, and therefore dif- 
ferential rotation generated in the collapse could persist. 
This "one-second-window" would provide sufficient time 
for several instabilities to develop in the PNS. Such insta- 
bilities are promising sources of gravitational waves. 

The restriction to the passive magnetic field approx- 
imation in studying magneto-rotational core collapse of 
weakly magnetized progenitors can be justified if the MRI 
is indeed inefficient, since none of the other estimated 
mechanisms seem to be able to amplify the magnetic 
field significantly on dynamical time scales. Otherwise, 
an "active" magnetic field approach becomes necessary. 
However, it has to be stressed that the use of active mag- 
netic fields alone for core collapse simulations will prob- 
ably not be sufficient to model all the effects amplifying 
the magnetic field, since the numerical resolution needed 
to correctly describe them (probably less than 10 m) is 
not affordable in current numerical simulations. In addi- 
tion most of the prospectively relevant effects have to be 
investigated in three dimensions, which makes the compu- 
tational task even more challenging. 
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Appendix A: Code tests 

Here wc discuss several tests we have designed in order to 
check the accuracy of our numerical code when solving the 
induction equation with the numerical methods described 
in this paper (see also Cerda-Duran & Font 2007). The 
"toroidal test" is set up for assessing the ability of the code 
to maintain various magnetic field configurations in equi- 
librium (labelled TTA and TTB) and to correctly com- 
pute the amplification of the toroidal magnetic field as it 
is wound up by a rotating fluid (TTC). On the other hand, 
the "poloidal test" (PT) is designed to check whether 
the code can correctly compute the compression of the 
poloidal magnetic field in a spherical collapse. Finally, the 
strong spherical explosion test checks that the code is able 
to handle the presence of radial shocks. We refer the in- 
terested reader to Cerda-Duran & Font (2007) for details 
on the setup of the toroidal and poloidal tests as well as 
the diagnostics we use to compute the errors and order of 
convergence of our numerical schemes. 

A.l. Toroidal tests 

Fig. A.l shows the global error a in the toroidal mag- 
netic field B v for the three tests TTA, TTB, and TTC 
against 1/f and the corresponding fits to a power law. 
Here / denotes the factor which specifies the increase in 
resolution from a coarse reference grid (see Cerda-Duran 
& Font 2007, for details). The resulting convergence or- 
der of each numerical scheme (minmod, MC, and PHM) 
as well as the errors for the highest resolution grid can 
be found in Table A.l. Our results show that (i) the or- 
der of convergence and the error is almost independent 
of the cell reconstruction scheme employed, (ii) the order 
of convergence for the TTC test is smaller than for the 
TTA and TTB tests, and (iii) the order of convergence for 
the tests TTA and TTB is N > 2, and hence higher than 
the theoretical expectation (which is second order, since 
it is limited by the order of the time discretization, for 
which we use a conservative, second order Runge-Kutta 
scheme) . 

The numbers reported in Table A.l demonstate that 
we obtain similar results in all three tests for linear recon- 
struction schemes (minmod and MC) and for the third 
order reconstruction scheme (PHM), as the order of the 
scheme is limited by the second order discretization in 
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1// 

Fig. A.l. Global error a in the toroidal magnetic field B v after 
a time evolution of 1 ms as a function of 1// for a sequence of 
models with grid resolutions 80 x 10 (/ = 1), 160 x 20 (/ = 2) 
and 320 x 40 (/ = 4). The top panel shows the error for the 
TTC test using different reconstruction schemes and the corre- 
sponding best fits to a power law: minmod (open circles, dotted 
line), MC (filled squares, dashed line), and PHM (filled circles, 
solid line). The bottom panel shows the error and respective 
fits using the PHM reconstruction for the TTA test (open cir- 
cles, dotted line), TTB (filled circles, dashed line), and TTC 
(open squares, solid line). 

Table A.l. Convergence order N for the tests performed 
(TTA, TTB, TTC, and PT) and for different reconstruction 
procedures (minmod, MC, and PHM). The error 0-320x40 for 
the higher resolution grid is also given. 



Test 


Reconstruction scheme 


N 


0320x40 


TTA 


minmod 


2.45 


1.2 x 10~ b 


TTA 


MC 


2.16 


2.4 x 10" 6 


TTA 


PHM 


2.46 


2.1 x 10 -6 


TTB 


minmod 


2.64 


7.7 x 10 _b 


TTB 


MC 


2.53 


1.2 x 10~ 5 


TTC 


PHM 


2.85 


0.5 x 10" 6 


TTC 


minmod 


1.38 


4.0 x 10" 5 


TTC 


MC 


1.48 


7.0 x 10 -5 


TTC 


PHM 


1.39 


3.5 x 10 -5 


PT 


minmod 


1.41 


8.3 x 10" 4 


PT 


MC 


1.11 


8.6 x 10" 4 


PT 


PHM 


1.17 


8.6 x 10" 4 



time and by the linear interpolation of the cell-centered 
magnetic fluxes (which is a consequence of using a stag- 
gered grid in the flux-CT scheme for the magnetic field). 



4r I I 1.0 ■ 

4 , , , °o 5 o ■ 

10 20 30 40 50 60 70 

Fig. A. 2. Local order of convergence (color coded) for the TTA 
test after a total time evolution of 1 ms. White color is used 
for values > 3.0. The horizontal and vertical axes represent the 
number of cells of the reference grid in the radial and angular 
direction, respectively. 

10* | | 
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t 

Fig. A. 3. Global error on r D* / B* 9 for the poloidal test (PT) 
as a function of time for three different grid resolutions: 80 x 10, 
160 x 20, and 320 x 40, and for four different reconstruction 
schemes (minmod, MC, and PHM). 



To understand these results we note that in test TTC there 
is a component of the magnetic field, B* v , which grows 
linearly in time, while in tests TTA and TTB no compo- 
nents evolve. Hence, the order of convergence for the latter 
is higher than for test TTC. This can be explained by in- 
vestigating the local order of convergence, i.e. the order 
obtained when computing the error a%j in each numerical 
cell instead of the global error a. The results for test TTA 
are displayed in Fig. A. 2 (similar plots can be obtained 
for the other two cases). At some particular grid zones 
the order of convergence is larger than two, while at most 
locations it remains around two. 

A.2. Poloidal test 

As mentioned before the setup and specifications of the 
poloidal test are described in detail in Cerda-Duran & 
Font (2007). Here we simply focus on showing the com- 
parison and performance of the various numerical schemes 
employed in our simulations. (Note that in Cerda-Duran 
& Font 2007, only the minmod scheme was assessed.) 
Fig. A. 3 shows the evolution of the error in the quan- 
tity r D* / B* at the equatorial plane (which is a quan- 
tity that should not change with time with respect to a 
Lagrangian coordinate system) during the spherical col- 
lapse of a 4/3-polytrope for different {r, 9} grid resolutions 
(80 x 10, 160 x 20, and 320 x 40), equally-spaced in the 
angular direction and logarithmically spaced in the radial 
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Fig. A. 4. Spherical explosion test at t = 4. The Lorentz factor 
W is color coded, and magnetic field lines are overplotted. 

direction. Table A.l gives again numbers for the error and 
the order of convergence of the various schemes computed 
at the end of the simulation (t = 20 ms). In all cases the 
errors are below 1%, even for the coarsest grid, and the 
order of convergence is higher than 1 (the presence of lo- 
cal extrema in the radial profiles of some hydrodynamical 
variables explains the reduction of the theoretical order 
as a built-in feature of total- variation diminishing numer- 
ical schemes) . Comparisons between the HLL approximate 
Riemann solver and the KT symmetric scheme yield al- 
most identical results (in agreement with Lucas-Serrano 
et al. 2004; Shibata & Font 2005; Anton et al. 2006). 

A. 3. Strong spherical explosion 

Explosions are among the most demanding tests for multi- 
dimensional codes as they show the ability of numerical 
schemes to handle shocks. Since the majority of existing 
MHD codes are written in Cartesian coordinates, the most 
common test is the cylindrical explosion. For relativistic 
MHD codes the setup of Komissarov (1999) for this test 
has been used by other authors (Del Zanna et al. 2003; 
Leismann et al. 2005) to compare different codes. However, 
in spherical coordinates it is not possible to impose the 
symmetries needed for this test. The most natural choice 
is thus the spherical explosion. Kossl et al. (1990) per- 
formed this test in the case of Newtonian MHD. To our 
knowledge, no spherical explosions tests have been per- 
formed in relativistic MHD. Therefore, we have designed 
such a spherical explosion test in which the initial jump 
conditions in the variables are the same as for the test 
by Komissarov (1999). In this way a relativistic shock is 
formed which does not occur in the Newtonian case of 
Kossl et al. (1990). 

Our test setup consists of an initial explosion zone with 
P = 1.0 and p = 10 -2 for r < 0.8, surrounded by an am- 
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Fig. A. 5. Results for the spherical test explosion at t = 4. The 
plots show profiles for the fluid pressure P (top panel) and the 
magnetic field component B° (bottom panel) along the equa- 
torial plane using different reconstruction schemes: minmod 
(dotted line), MC (dashed line), and PHM (solid line). The 
grid resolution is 320 x 40. 

bicnt gas with P = 3 x 10~ 5 and p = 10~ 4 . The explosion 
region joins the ambient medium by matching an expo- 
nential decline in a transition region region 0.8 < r < 1.0. 
The velocities are initially zero, and the magnetic field 
is homogeneous and parallel to the symmetry axis. The 
background spacetime is considered to be flat. The inital 
data are evolved using an ideal gas EOS with adiabatic 
index 7 = 4/3. We use an evenly spaced grid with a max- 
imum radius of r = 6.0. We perform the test for three 
resolutions (80 x 10, 160 x 20, and 320 x 40) for all recon- 
struction schemes. 

Fig. A. 4 shows the Lorentz factor W at t — 4.0. A 
strong spherical shock has formed, propagating close to 
the speed of light, and as a consequence the magnetic field 
lines are compressed in the direction perpendicular to the 
axis. The results for this test are qualitatively compara- 
ble to the weakly magnetized case in Komissarov (1999). 
Fig. A. 5 shows radial profiles for P and B 6 along the equa- 
torial plane at the end of the simulation, using various 
reconstruction schemes. These plots are qualitatively sim- 
ilar to those of the cylindrical explossion (see e.g. Fig. B.4 
in Leismann et al. 2005). All numerical schemes exhibit 
first order convergence with increasing resolution, as is ex- 
pected to happen at shocks. The MC and PHM schemes 
yield very similar results, while minmod shows slightly 
lower values. 
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Appendix B: Estimation of the growth rates of 
the fi-dynamo 

To compute the characteristic time scales on which the 0- 
dynamo mechanism amplifies the magnetic field one has 
to study how the wind-up proceeds. Let us consider a sta- 
tionary rotating configuration with no meridional flows, 
w'^^^Oandw*^ Q*(r,6)r sin (9, where fi*(r,0) 
stands for the rotation law. Under these conditions and 
in the passive field approximation, the induction equation 
can be integrated analytically. The solution shows that the 
poloidal component of the magnetic field remains constant 
and the toroidal component grows linearly with time as 

B* v (t) = B* v (t = 0) + tmB* ■ Vn*. (B.l) 

This equation specifies the toroidal magnetic field at any 
given time, provided that the poloidal component is con- 
stant and the angular velocity profile is fixed. For a time 
t > B'^ a /(zuB* ■ VST), which is ~ 1 ms in our simula- 
tions, we can use this expression to compute the magnetic 
energy 

E magv « J d 3 x^- (Tu|Vfi*|) 2 t\ (B.2) 

where is the component of B* parallel to Vfl*. The 
rotation profiles of the final PNS can be approximated 
in all our models by the rotation law (25) (Villain ct al. 
2004). In the Newtonian limit (26), which is good enough 
for this estimate, we can compute an upper limit to 
the magnetic energy considering the maximum value of 
|o7VO*| max = fi*/2, which yields 

E maglp < BmagO t ■ (B.3) 

Therefore, an estimate for the upper limit of the amplifi- 
cation of the magnetic energy parameter is 

n* 2 n n* 2 „ / t \ 2 

/3 mag « V < Pn^-t 2 < /3 P oio^ t 2 = I- , (B.4) 

where we have defined the time scale for amplification of 
the magnetic field by the Jl-dynamo as 

m = 1 (B.5) 

V Ppolo 

This gives us the characteristic time scale in which /3 mag 
reaches a value of 1; therefore, the saturation time i sat 
should be a fraction of this time. As the Sl-dynamo oper- 
ates by transforming rotational energy into magnetic en- 
ergy, the maximum energy can be extracted by the mag- 
netic field is the one that is contained in the differential 
rotation of the core. In accordance with numerical simu- 
lations using strong magnetic fields (Obergaulinger et al. 
2006b) we estimate this amount to be 10% of the total 
rotational energy, i.e. (3 mag (t sat ) = 0-1 (3 Iot (t sat ). We also 
assume that the evolution of the magnetic energy parame- 
ter is given by Eq. (B.4) and that the energy is conserved, 

i.e. A-ot(i) = /?rot(*0) - AnagW- 



